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(57) Abstract: Blind source separation (BSS) of statistically independent signals with low signal-to-noise plus interference ratios 



under a narrowband assumption utilizes cumulants in conjunction with spectral estimation of the signal subspace to perform the blind 
separation. The BSS technique utilizes a higher-order statistical method, specifically fourth-order cumulants, with the generalized 
eigen analysis of a matrixpencil (13) to blindly separate a linear mixture of unknown, statistically independent, stationary narrow- 
band signals at a low signal-to-noise plus interference ratio having the capability to separate signals in spatially and/or temporally 
correlated Gaussian noise. This BSS provides the ability to blindly separate signals in situations where no second-order technique 
has been found to perform the blind separation satisfactorily, for example, at a low signal-to-noise ratio when the number of sources 
equals the number of sensors or when the noise is spatially and temporally colored. 
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BLIND SOURCE SEPARATION UTILIZING A SPATIAL FOURTH ORDER 

CUMULANT MATRIX PENCIL 

The present application claims priority to U.S. Provisional Application Serial 
No. 60/374,149 filed April 22, 2002, entitled "Blind Source Separation Using A 
5 Spatial Fourth Order Cumulant Matrix Pencil" the entirety of which is incorporated 
herein by reference. 

The present invention is generally related to separating individual source 
signals from a mixture of source signals, and more specifically related to blind source 
separation. 

0 A classic problem in signal processing, often referred to as blind source 

separation (BSS), involves recovering individual source signals from a composite 
signal comprising a mixture of those individual signals. An example is the familiar 
"cocktail party" effect, wherein a person at a party is able to separate a single voice 
from the combination of all voices in the room. The separation is referred to as 

5 "blind" because it is often performed with limited information about the signals and 
the sources of the signals. 

Blind source separation (BSS) is particularly applicable to cellular and 
personal wireless communications technologies, wherein many frequency bands have 
become cluttered with numerous electromagnetic emitters, often co-existing in the 

10 same spectrum. The problem of co-channel emitters is expected to only worsen in 
years to come with the development of low power, unlicensed wireless technologies 
such as Bluetooth® and other personal area networks. These developments have 
resulted in the use of multiple sensors and array signal processing techniques to 
perform spectral monitoring. Such techniques enable the exploitation of spatial 

[5 information to separate co-channel emitters for detection, classification, and 

identification. Additionally, many signals designed for a low probability of detection 
(LPD) or low probability of intercept (LPI) may use ambient background 
electromagnetic radiation and known co-channel emitters as a means of concealment. 
Constructing single sensor receiver systems with the required sensitivity to such 

\0 emitters is generally prohibitive. Thus, many applications utilize BSS and sensor 
arrays. 

Several techniques have been proposed to solve the BSS problem. These can 
be classified into two main groups. The first being based on second-order statistics, 
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and the second being based on higher-order statistics, such as those based on 
independent components analysis (ICA) and other higher-order spectral estimation and 
spatial filtering techniques. 

One second-order blind source separation technique is a spectral estimation 
5 method that exploits the rotational invariance of the signal subspace to estimate the 
direction of arrival. This technique known as Estimation of Signal Parameters via 
Rotational Invariance (ESPRIT) employs pairs of calibrated elements and uses a 
matrix-pencil formed by the spatial correlation and cross-correlation matrices. See, for 
example, R Roy, A. Paulraj, T. Kailath, "Direction-of- Arrival Estimation by Subspace 

10 Rotation Methods,'* Proa ICASSP86, pp. 2495-2498 and R. Roy and T. Kailath, 
"ESPRIT - Estimation of Signal Parameters via Rotational Invariance Techniques," 
IEEE Trans, on ASSP, Vol. 37, No. 7, Jul. 1989, pp. 984-995, which are each 
incorporated by reference in their entirety as if presented herein. However, a 
disadvantage of ESPRIT is that at very low signal-to-noise ratios the signal plus noise 

15 subspace and the noise subspace are indistinguishable thus making the conversion of 
the noise subspace of the spatial correlation matrix impractical. This is due in part to 
ESPRIT requiring an estimation of the noise variance to convert the noise subspace 
into a null subspace of the spatial correlation matrix, and assuming that the noise is 
spatially white. 

20 Another second order blind source separation technique is known as the 

Constant Modulus Algorithm (CMA), also referred to as Goddard's algorithm. The 
CMA is an adaptive spatial filtering technique, which is used to perform source 
separation by determining a set of spatial filter tap weights that forces the output signal 
to have a modulus as close to unity as possible. Typically, the CMA is performed 

25 sequentially to separate all source signals. The CMA has been suggested for use as a 
blind equalization technique to reduce inter-symbol interference of constant modulus 
signals, such as FSK, PSK, and FM modulated signals, on telephone channels (See, for 
example, D.N. Godard, "Self-recovering Equalization and Carrier Tracking in Two- 
dimensional Data Communication Systems," IEEE Trans. Commun., Vol. COMM-28, 

30 Nov. 1980, pp. 1867-1875, which is incorporated by reference in its entirety as if 

presented herein.), and to perform blind equalization to combat multi-path fading and 
to suppress co-channel interfering signals (See, for example, B.G. Agee, "The 
Property Restoral Approach to Blind Adaptive Signal Extraction," Ph.D. Dissertation, 
Dept. Elect. Eng. Comput. Science, Univ. of Calif, Davis, 1989, which is incorporated 
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by reference in its entirety as if presented herein). However, the CMA technique 
works only for signals with a constant modulus and is not practicable for most 
applications. In practice, because signals are filtered to limit their spectral occupancy 
at the transmitter and to limit the noise bandwidth at the receiver, true constant 
5 modulus signals rarely exist. Furthermore, at very low signal-to-noise ratios, noise 
dominates the input signal thus distorting the spatial filter's output signal's modulus 
and causing large fluctuation in the error signal used in adaptation. 

Yet another second order blind source separation technique is a spatial filtering 
technique using second-order cyclostationary statistics with the assumption that the 

10 source signals are cyclostationary. This technique was developed as a blind single- 
input single-output (SISO) channel identification technique for use in blind 
equalization (See, for example, L. Tong, G. Xu, and T. Kailath, "Blind Identification 
and Equalization Based on Second-Order Statistics: A Time-Domain Approach," 
IEEE Trans. Information Theory, Vol. 40, No. 2, Mar. 1994, pp. 340-349, which is 

1 5 incorporated by reference in its entirety as if presented herein), and was later adapted 
to perform the blind separation of cyclostationary signals (See, for example, L. 
Castedo and A.R. Figueiras-Vidal, "An Adaptive Beamforming Technique Based on 
Cyclostationary Signal Properties," IEEE Trans. Signal Processing, Vol. 43, No. 7, 
Jul. 1995, pp. 1637-1650, which is incorporated by reference in its entirety as if 

20 presented herein). One disadvantage of this cyclostationary approach is that it requires 
different symbol rates and/or different carrier frequencies for seperating multiple 
superimposed signals. Another disadvantage is that residual carrier offsets with a 
random initial phase can cause the signals to become stationary, causing the 
cyclostationary assumption to become invalid. Other disadvantages include the fact 

25 that this approach precludes separating sources that may use non-linear or non-digital 
modulations, and this approach assusmes the noise vector is temporally and spatially 
white. 

Still another blind source separation technique based on second-order statistics 
is referred to as Second-Order Blind Identification (SOBI). See, for example, A. 
30 Belouchrani, K. Abed-Meraim, J.F. Cardoso, and E. Moulines, "Blind Source 

Separation Using Second-Order Statistics," IEEE Trans. Signal Processing, Vol. 45, 
No. 2, Feb. 1997, pp. 434-444, for a description of this technique, which is 
incorporated by reference in its entirety as if presented herein. This technique exploits 
the time coherence of the source signals and relies on the joint diagonalization of a set 
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of covariance matrices. A disadvantage of this technique is that it requires additive 
noise to be temporally white and uses the eigenvalues of the zero lag matrix to 
estimate the noise variance and to spatially whiten the sensor output vector. Another 
disadvantage is that at low signal to noise rations, the estimation of the noise variance 
is difficult at best and impossible in most cases. Yet another disadvantage is that the 
number of sources must be known or estimated. Still, another disadvantage is that the 
SOBI technique is valid for spatially correlated noise. The estimation of the noise 
covariance is extremely difficult even at high signal-to-noise ratios, thus making the 
technique impracticable. 

Another second-order blind source separation technique is based on the 
generalized eigen decomposition of a spatial covariance matrix-pencil. This technique 
is related to the ESPRIT algorithm, but does not require 2N sensors to separate up to 
N-l signals because time and/or polarization diversity is used to estimate the pair of 
spatial covariance matrices. This technique uses a dual polarized array and estimates 
the spatial covariance matrices on the two orthogonal polarizations to form the matrix 
pencil. See, for example, A. Belouchrani, K. Abed-Meraim, J.F. Cardoso, and E. 
Moulines, "Blind Source Separation Using Second-Order Statistics," IEEE Trans. 
Signal Processing, Vol. 45, No. 2, Feb. 1997, pp. 434-444, which is incorporated by 
reference in its entirety as if presented herein. This later evolved into using spatial 
covariance matrices at a zero time lag and a non-zero time lag to form the matrix- 
pencil. One disadvantage of this technique is that it is limited to separating up to N-l 
sources with N sensors. This is due in part to the approach requiring the estimation of 
the noise variance, similar to ESPRIT, and assuming that the noise is spatially and 
temporally white. 

Finally, another second order blind source separation technique utilizes two 
non-zero time lags in the estimation of the spatial covariance matrices. See, for 
example, C. Chang, Z. Ding, S. F. Yau, and F. H. Y. Chan, "A Matrix-Pencil 
Approach to Blind Separation of Non- White Sources in White Noise," Proc. 
ICASSP98, Vol. IV, pp. 2485-2488, and C. Chang, Z. Ding, S. F. Yau, and F. H. Y. 
Chan, "A Matrix-Pencil Approach to Blind Separation of Colored Non-Stationary 
Signals/' IEEE Trans. Signal Processing, Vol. 48, No. 3, Mar. 2000, pp. 900-907, 
which are each incorporated by reference in their entirety as if presented herein. The 
non-zero time lags combined with the assumption that the noise vector is temporally 
white eliminates the need to estimate the noise variance(s) in order to remove the noise 
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subspace and thus allows the technique to separate up to N sources with N sensors. 
However, disadvantages of the this second-order matrix-pencil technique include the 
requirement that the noise vector be temporally white and the fact that the estimate of 
the separation matrix is not done with one of the spatial covariance matrices at a zero 
5 time lag, which is when the signal auto-correlations are at their maximum values. 

These disadvantages are exacerbated by the fact that in many practical applications the 
noise bandwidth is limited to be on the order of the signal bandwidth making the noise 
and signal decorrelation times to be approximately equal. 

Higher-order blind source separation techniques include all methods that 
10 employ statistics of order greater than two. These include independent component 
analysis (ICA) methods, spatial filtering methods, and spectral estimation based 
methods and can use either moments or cumulants of order three or higher. 

The independent components analysis (ICA) methods seek a separating matrix 
that maximizes the statistical independence of the outputs of the separation process. 

1 5 See, for example, A.K. Nandi, Blind Estimation Using Higher-Order Statistics. 
(Kluwer Academic, Dordecht, The Netherlands: 1999), A. Hyvarinen, "Survey on 
Independent Component Analysis," Neural Computing Surveys, Vol. 2, No. 1, 1999, 
pp. 94 - 128, and J.F. Cardoso, "Blind Signal Separation: Statistical Principles"' Proc. 
of the IEEE, Vol. 9, No. 10, Oct. 1998, pp.2009-2025, which are each incorporated by 

20 reference in their entirety as if presented herein. Various ICA methods are known, and 
are primarily differentiated by the objective/contrast function used to measure 
statistical independence. Examples of ICA based blind source separation algorithms 
include the Jutten-Herault algorithm (See, for example, C. Jiitten and J. Herault, 
"Blind Separation of Sources, Part I: An Adaptive Algorithm Based on Neuromimetic 

25 Architectures," Signal Processing, Vol. 24, 1991, pp. 1-10, which is incorporated by 
reference in its entirety as if presented herein), which attempts to achieve separation 
via canceling non-linear correlations through the use of a neural network; the higher- 
order eigenvalue decomposition or HOEVD method (See, for example, P. Comon, 
"Independent Component Analysis, A New Concept?," Signal Processing, Vol. 36, 

30 No. 3, Apr. 1 994, pp. 287-3 14, which is incorporated by reference in its entirety as if 
presented herein); the joint approximate diagonalization of eigen matrices (JADE) 
algorithm (See, for example, J.F. Cardoso and A. Souloumiac, "Blind Beamforming 
for Non-Gaussian Signals " IEE Proceedings F, Vol. 140, No. 6, Dec. 1993, pp. 362- 
370, which is incorporated by reference in its entirety as if presented herein), which 
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exploits the eigen structure of the fourth-order cumulant tensor; the information 
maximization or infomax technique (See, for example, A J. Bell and T.J. Sejnowski, 
"An Information-Maximization Approach to Blind Source Separation and Blind 
Deconvolution " Neural Computing, Vol. 7, 1995, pp. 1 129-1 159, which is 
5 incorporated by reference in its entirety as if presented herein), which seeks a 

separation matrix that maximizes the output entropy; and the equivariance adaptive 
source separation or EASI algorithm (See, for example, J.F. Cardoso and B. Hvam 
Laheld, "Equivariance Adaptive Source Separation," IEEE Trans. On Signal 
Processing, Vol. 44, No. 12, Dec. 1996, pp. 3017-3030, which is incorporated by 

10 reference in its entirety as if presented herein), in which an estimate of the mixing 

matrix is chosen such that transformations in the sensor output data produces a similar 
transformation in the estimated mixing matrix. Disadvantages of these ICA 
techniques include (1) a pre- whitening step is required, (s), the techniques need to be 
parametrically timed to the application, and (3) the ICA based techniques tend to have 

1 5 a slow convergence time. 

Another higher-order blind source separation technique, known as the Kurtosis 
Maximization Algorithm (KMA), utilizes kurtosis as a separation measure (See, for 
example, Z. Ding, "A New Algorithm for Automatic Beamforming," Proc. 25 th 
Asilomar Conf. Signals, Syst., Comput, Vol. 2, 1991, pp. 689-693, and Z. Ding and T. 

20 Nguyen, "Stationary Points of Kurtosis Maximization Algorithm for Blind Signal 

Separation and Antenna Beamforming," IEEE Trans. Signal Processing, Vol. 48, No. 
6, Jun. 2000, pp. 1587-1596, which is incorporated by reference in its entirety as if 
presented herein). The KMA is an adaptive spatial filtering technique with a non- 
deterministic convergence time. One disadvantage of the KMA is that it can not 

25 simultaneously separate the source signals. The KMA requires that the sources be 
separated sequentially one at a time, starting with the signal having the largest 
kurtosis. Thus, no knowledge of the number of signals is provided by the technique. 
Other disadvantage of the KMA are that it requires the noise to be spatially white and 
its convergence with multiple signals has only been proven in the noise free case. 

30 Finally, another higher-order BSS technique is a higher-order version of the 

ESPRIT algorithm. As the name implies, higher-order ESPRIT replaces the spatial 
correlation or covariance matrices with spatial fourth-order cumulant matrices. See 
H.H. Chiang and C.L. Nikias, "The ESPRIT Algorithm with Higher-Order Statistics," 
Proc. Workshop on Higher-Order Spectral Analysis, Vail, CO., Jun. 1989, pp. 163- 
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168, C.L. Nikias and A.P. Petropulu, Higher-Order Spectra Analysis: A Non-Linear 
Signal Processing Framework. (PTR Prentice-Hall, Upper Saddle River, NJ: 1993), 
and M.C. Dogan and J.M. Mendel, "Applications of Cumulants to Array Processing - 
Part I: Aperture Extension and Array Calibration," IEEE Trans. Signal Processing, 
5 Vol. 43, No. 5, May 1995, pp. 1200-1216, for descriptions of various types of higher- 
order ESPRIT techniques, which are each incorporated by reference in their entirety as 
if presented herein. These higher-order ESPRIT techniques possess several 
disadvantages. Some of these higher-order ESPRIT techniques require the calibration 
of the N pairs of sensors but can now separate up to N sources since the noise variance 
10 no longer needs to be estimated. Other higher-order ESPRIT techniques require the 
array to be calibrated to guarantee that the pairs of sensors have identical manifolds 
(similar to the standard ESPRIT). These techniques degrade in performance 
significantly as the sensor pairs' manifolds deviate from one another. 

Each of the above-mentioned blind source separation techniques have the 
15 disadvantages noted. Additionally, none of the above-mentioned blind source 

separation techniques operate satisfactorily in the situation where there is a low signal- 
to-noise plus interference ratio. Accordingly, an improved blind source separation 
technique is desired. 

In one embodiment of the present invention, a method for separating a plurality 
20 of signals provided by a respective plurality of sources and received by an array 
comprising a plurality of elements, includes generating a separation matrix as a 
function of time differences between receipt of the plurality of signals by the plurality 
of elements and a spatial fourth order cumulant matrix pencil. The method also 
includes multiplying the separation matrix by a matrix representation of the plurality 
25 of signals. 

In another embodiment of the present invention, a system for separating a 
plurality of signals provided by a respective plurality of sources includes a receiver for 
receiving the plurality of signals and for providing received signals. The system also 
includes a signal processor for receiving the received signals, generating a separation 
30 matrix, and multiplying the separation matrix by a matrix representation of the 
received signals. The separation matrix is a function of time differences between 
receipt of the plurality of signals by the receiver and a function of a spatial fourth 
order cumulant matrix pencil. 
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Brief Description Of The Drawings 

In the drawings: 

Figure 1 is a functional block diagram of a system for performing blind source 
separation utilizing a spatial fourth order cumulant matrix pencil in accordance with an 
5 embodiment of the present invention; 

Figure 2 is an illustration of signal source, array elements, and a processor for 
performing array signal processing and BSS processing in accordance with an 
embodiment of the present invention; 

Figure 3 is an illustration of a MIMO blind channel estimation scenario 
1 0 showing five unknown sources having distinct radiating patterns and five sensors 
having distinct receiving patterns; 

Figure 4 is a graphical illustration of time delays between sensors and sources; 

Figure 5 is an illustration depicting blind source separation (BSS) showing an 
input signal mixed with noise provided to the separation process; 

1 5 Figure 6 is an illustration depicting repeating the separation process for a single 

repeated eigenvalue; 

Figure 7 is a flow diagram of a process for performing blind source separation 
using the spatial fourth-order cumulant matrix-pencil in accordance with an 
embodiment of the present invention; and 

20 Figure 8 is a continuation of the flow diagram of Figure 7. 

Detailed Description 

A technique for performing blind source separation (BSS) in accordance with 
the present invention utilizes cumulants in conjunction with spectral estimation of the 

25 signal subspace to perform the blind separation of statistically independent signals 
with low signal-to-noise ratios under a narrowband assumption. This BSS technique 
makes use of the generalized eigen analysis of a matrix-pencil defined on two similar 
spatial fourth-order cumulant matrices. The herein described BSS technique utilizes a 
higher-order statistical method, specifically fourth-order cumulants, with the 

30 generalized eigen analysis of a matrix-pencil to blindly separate a linear mixture of 
unknown, statistically independent, stationary narrowband signals at a low signal-to- 
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noise ratio having the capability to separate signals in spatially and/or temporally 
correlated Gaussian noise. This BSS technique provides a method to blindly separate 
signals in situations where no second-order technique has been found to perform the 
blind separation, for example, at a low signal-to-noise ratio when the number of 
5 sources equals the number of sensors. 

To describe this BSS technique, a definition of a spatial fourth-order cumulant 
matrix suited to blind source separation with non-equal gain and/or directional sensors 
and a definition of a spatial fourth-order cumulant matrix-pencil using temporal 
information are provided. The herein description also utilizes the concept of 
10 separation power efficiency (SPE) as a measure of the BSS technique's performance, 
and applies the concept of wide sense equivalence between matrix-pencils to the field 
of matrix algebra. 

As an overview, the BSS technique described herein utilizes cumulants in 
conjunction with a spectral estimation technique of the signal subspace to perform 

1 5 blind source separation in the presence of spatially and/or temporally correlated noise 
at low signal-to-noise ratios. Prior to deriving a separation algorithm based on 
cumulants, a narrowband array model is developed, all assumptions are stated, four 
performance measures are defined, and the relevant cumulant properties that allow for 
the spatial mixing matrix information to be extracted from a spatial cumulant matrix 

20 are presented. A novel spatial cumulant matrix definition is then developed and its' 
relevant matrix properties are derived in order to determine which mathematical 
methods are valid for extracting the spatial information about the mixing matrix. 
Additionally, two alternative definitions for the spatial fourth-order cumulant matrix 
are described and relevant properties are derived. Furthermore, the definitions, 

25 properties, and use of a generalized eigen analysis of a matrix-pencil defined on two 
similar spatial fourth-order cumulant matrices are explored and their applicability to 
solving the blind source separation problem is investigated. A process is described for 
performing the blind source separation based on the signal subspace technique using 
matrix-pencils. In the process the concept of wide sense equivalence between matrix- 

30 pencils is developed and then used to show that the generalized eigenvalues of a 

matrix-pencil defined on two similar spatial fourth-order cumulant matrices are equal 
to the ratio of fourth-order cumulant of each source at a set of time lags (0,0,0) to the 

fourth-order cumulant at the set of lags, (r, , r 2 , r 3 ) . Thus the concept of a normalized 



9 



WO 03/090127 



PCTAJS03/11473 



fourth-order auto-cumulant function is introduced. To further aid in understanding 
this BSS technique, notation used herein is presented below. 

Ms Number of Sources 
N= Number of Sensors 

5 Pj s Normalized Power of the j th source signal 

ntj (/) a Continous Time Unit Power Modulated Signal from the 7 th source 

Sj (/) = Continous Time Signal from the 7 th source = (t) 

r j {t) = Delayed version of Sj (t) 

x, (t) = Continous Time Signal from the 1* sensor. 

x(/)s The vector of sensor outputs. 

h fJ (t) s Continous Time Impulse Response of the channel between the 7 * source 
and the 1* sensor 

n % (t) s Additive Noise Process at the i* sensor. 

a/ = Variance of the Noise Process at the i th sensor. 

r~ = Propogation Delay from the j mth source to the z lh sensor 

Ar, k j s "Differential Time Delay". The difference in propogation delay from the output 

of the 7 th source to the sensor output and from the output of the j th source 
to the /* sensor output. 

Tj = "Reference Time Delay" from the /^source to some arbitrary array reference point 
10 in the vicinity of the array. Nominally this can be the average propogation delay to 

all N sensors from the /* source. 

Azy = "Relative Time Delay". The difference in propogation time from the f h source 
to the i th sensor and the array refence point. 

x s Time Difference in Correlation of Stationary Processes 

v r = Complex Weight for the 7 th source at the i ,h sensor for the Narrow Band Model. 
tJ The ij element of the f, Mixing Matrix". The i element of the / steering vector. 

15 

Vj = The f "Steering Vector" for the Narrow Band Model. 
V = The Narrow Band Model "Mixing Matrix". 
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w fJ = Complex Weight for the j mfh source at the sensor for the Narrow Band Case. 

The ij element of the "Separating Matrix". The/ th element of the/ 11 sensor 
weight vector. 

W = The "Separation Matrix". 

a fj s Real valued gain(attenuation) of the channel from the i source output to the j 
sensor output. 

BW NEq [ ]s Noise Equivalent Bandwidth 
5 BW COH = Coherence bandwidth of the Channel between the /^source and the i ih sensor. 



y. (/) = The f h output from the separation process. It is a noisy estimate of the 
of the /^delayed source signal, r } (t). 

y (/) s The vector of output signals from the separation process. 

Pj s The /* signal loss term. Element of the "loss" matrix. 
10 S , s The separation process output signal power of the /* source signal. 

I j s The residual interference power in the j lh separation process output. 

Nj = The noise power in the f separation process output. 

£j s The "Interference -to- Signal Ratio " for the f separation process output. 

ISR^ s The "Average Interference - to - Signal Ratio". 
1 5 ISRun 55 The "Maximum Interference - to - Signal Ratio". 

£ . s The "Power Efficiency" of a blind source separation algorithm for the /* source. 

Jj = The "Average Power Efficiency" of a blind source separation algorithm. 

^min = Tne "Minimum Power Efficiency" of a blind source separation algorithm. 

C*(r„7- 2 ,r 3 )3 NxN "Spatial Fourth -Order Cumulant Matrix 1" with delay 

lags r„r 2 ,r 3 . 

20 C * ( r i ' T 2 » r 3 ) = NxN " Spatial Fourth - Order Cumulant Matrix 2 " with delay 

lags T t ,T 2 ,T y 
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C* (r,,r 2 ,r 3 )= NxN "Spatial Fourth-Order Cumulant Matrix 3" with delay 

lags r 15 r 2 ,r 3 . 

Cum[ ]s Cumulant Operator. 

c tj ( r i> r 2» r 3) - The fourth-order cumulant of the f h source signal with delay 

lags t {> t 29 t 3 . Also referred to as the fourth-order auto-cumulant. 

V s The "Modified Mixing Matrix", Defined as the Hadamard Product V © V © V. 

ZT r 4 (r p r 25 r 3 ) = The normalized fourth-order cumulant of the j lh source signal with 

5 delay lags r, , r 2 , r 3 . Also referred to as the normalized fourth-order 

auto-cumulant. 

C* (r p r 2 ,r 3 ) s MxM Diagonal " Four th - Order Signal Cumulant Matrix" with delay 

lags r l5 r 2 ,r 3 . 

C ( ) = "Column Space" of a matrix. 

N r ( ) a The "JJ/g/A M/// 5/7ace" of a matrix. 
10 , ( ) = The "Left Null Space" of a matrix. 

1^ == N x N Identity Matrix. 

tr[ )= The "Trace" of a matrix. 
)= The "Span" of a sub-space. 

p( )= The "fowfc" of a matrix. 

15 t = Vector notation for the set of delay lags, {r, , r 2 , r 3 } . 

P x (/1,t) = The "Spatial Fourth - Orsfer Cumulant Matrix - Pencil" using a pair of 
Spatial Fourth-Order Cumulant Matrix l's. 

P x 55 The "Spatial Fourth - Order Cumulant Matrix - Pencil" using a pair of 
Spatial Fourth-Order Cumulant Matrix 2's. 

P*(A,t)s The "Spatial Fourth - Order Cumulant Matrix - P<?/?cz7" using a pair of 
Spatial Fourth-Order Cumulant Matrix 3's. 

P r (A,t) s The "Fourth - Oder S/gna/ Cumulant Matrix - Pencil" using a pair of 
Diagonal Fourth-Order Signal Cumulant Matrices. 
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A(A,B) = The "Spectrum" of the pencil defined on the matrices A and B. 
The set of generalized eigenvalues. 

i(A,B) s The "Finite Spectrum" of the pencil defined on the matrices A and B. 
The set of non-zero finite generalized eigenvalues. 

A, s The "j' h Eigenvalue" of the pencil defined on a pair of spatial fouth-order cumulant 
matrices. There are M such eigenvalues, counting multiplicities. X } takes on one 
of the K values of fx k . 

ju k s The "k th Distinct Eigenvalue" of the pencil defined on a pair of spatial fouth-order 
cumulant matrices. There are K such values that the set of AJs takes on. 

g k = The set of indeices, { j} , where X } = fi k . 

e 7 = The Nxl "j lh Eigenvector" of the pencil defined on a pair of spatial fouth-order 
cumulant matrices associated with the eigenvalue X y 

10 y = The "Normalization Factor" for the y ,A eigenvector. y, = 

I J 

rjl eom = The "Geometric" Multiplicity of an Eigenvalue. 

r]f % s The "Algebraic" Multiplicity of an Eigenvalue. 

7j k = The "Multiplicity" of an Eigenvalue when 7jf om =rj k g . 

15 Figure 1 is a functional block diagram of a system 100 for performing blind 

source separation utilizing a spatial fourth order cumulant matrix pencil in accordance 
with an embodiment of the present invention. System 100 comprises a receiver 1 1 and 
a signal processor 12. The receiver 1 1 receives signal s(t), which is indicative of a 
plurality of signals provided by a respective plurality of sources and provides signal 

20 x(t) to the signal processor 12. The receiver 1 1 may be any appropriate receive 

configured to receive the signal s(t). For example, the signal s(t) may be an acoustic 
signal, an optical signal, a seismic signal, an electromagnetic signal, or a combination 
thereof, and the receiver 1 1 may be configured to receive the respective type of signal. 

* 
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In one embodiment, the receiver 1 1 is configured as an array having a plurality of 
elements. The signal s(t) is received and appropriately processed (e.g., time delayed 
and multiplexed) and provided to the signal processor 14 in the form of signal x(t). 

The signal processor 12 may be any appropriate processor configured to 

5 processor the signal x(t), such a general purpose computer, a laptop computer, a 
special purpose computer, a hardware implemented processor, or a combination 
thereof. The signal x(t) may be in any appropriate format, such as an optical signal, 
and electromagnetic signal, a digital signal, and analog signal, or a combination 
thereof. As will be explained in more detail below, the signal processor 12 comprises 

10 a matrix pencil estimation portion 13, a non-zero finite eigenvalue determination 
portion 14, a number of distinct eigenvalues determination portion 15, a multiplicity 
determination portion 16, a linearly independent eigenvector calculation portion 17, a 
normalization factor calculation 1 8, a separation vector generation portion 1 9, a 
separation matrix generation portion 20, and an optional separation power efficiency 

1 5 calculation portion 2 1 . The matrix pencil estimation portion 1 3 is configured to 
estimate the spatial fourth order cumulant matrix pencil as a function of time 
differences of the arrival of the signal s(t) at the elements of the receiver 1 1 . The non- 
zero finite eigenvalue determination portion 14 is configured to determine the non- 
zero finite eigenvalues for the spatial fourth order cumulant matrix pencil. The 

20 number of distinct eigenvalues determination portion 15 is configured to determine the 
number of eigenvalues that are distinct. The multiplicity determination portion 16 is 
configured to determine the multiplicity of each of the distinct finite eigenvalues. The 
linearly independent eigenvector calculation portion 1 7 is configured to calculate 
linearly independent eigenvectors for each of the distinct finite eigenvalues. The 

25 normalization factor portion 1 8 is configured to calculate, for each eigenvalue having a 
multiplicity equal to one, a normalization factor and to generate a respective separation 
vector as a function of the normalization factor and an eigenvector corresponding to 
the eigenvalue having a multiplicity equal to one. The separation vector generation 
portion 19 is configured to generate, for each repeated eigenvalue, a separation vector 

30 as a function of an eigenvector corresponding to the repeated eigenvalue. The 

separation matrix generation portion 20 is configured to generate the separation matrix 
as a function of the separation vectors. The optional separation power efficiency 
calculation portion 21 is configured to calculate the efficiency of the separation 
process in accordance with the following formula: ^pS/Pj, wherein CJj is indicative of 
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the separation power efficiency for the j source of the plurality of sources, Sj is 
indicative of a power of a separated signal from the 3 th source, and Pj is indicative of a 
normalized power of a signal from the j 1 * 1 source. 

Figure 2 is an illustration of signal source 24, array elements 26, and a 

5 processor 22 for performing array signal processing and BSS processing in accordance 

with an embodiment of the present invention. Array signal processing is a 

specialization within signal processing concerned with the processing of a set of 

signals generated by an array of sensors at distinct spatial locations sampling 

propagating wavefields, such as electromagnetic, seismic, acoustic, optical, 

10 mechanical, thermal, or a combination thereof, for example. As shown in Figure 2, the 
array samples the f h wavefield, r j (f,z,) , generated by the f source 24j at locations 

{z x ,z 2 ,. .. >z N ) (only one location, Zj, shown in Figure 2) with a set of sensors 26 j 

which generate signals Xj(t) indicative of the wavefield at each location, Zj. The 

signals Xi(t) may be any appropriate type of signal capable of being processed by the 

1 5 processor 22. Examples of appropriated types of signals Xj(t) include electrical 

signals, acoustic signals, optical signals, mechanical signals, thermal signals, or a 
combination thereof. The signal x, (t) provided by the i th sensor, 26j, comprises the 

sum of the wavefields from all sources 24 at each sensor's location, each weighted 
with response of the sensor in the signal's r j ) direction of arrival, plus an 

20 additive noise term, /?, (t) . As described in more detail herein, the processor 22 

processes the signals x(t) for enhancing sets of sources signals' individual signal-to- 
interference-plus-noise ratios by suppressing interfering source signals at different 
spatial locations without knowledge of the source signal characteristics, the channels 
between the sources and the array elements, the sources' locations, or array geometry 
25 via a blind source separation (BSS) technique in accordance with the present 
invention. 

A blind source separation technique in accordance with the present invention is 
described herein by defining underlying assumptions made about the source signals 
and noise sources. Different multiple input multiple output (MIMO) array channel 
30 models are described resulting in a narrowband model, which is utilized in the BSS 
technique in accordance with present invention. 

Blind source separation (BSS) is applicable to many areas of array signal 
processing that require the enhancement and characterization of an unknown set of 
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source signals generated by a set of sensors that are each a linear mixture of the 
original signals. These include, for example, signal intelligence, spectral monitoring, 
jamming suppression, and interference rejection, location, and recognition. Typically, 
the mixing transformation, source signal characteristics, and sensor array manifold are 
5 unknown. Thus, blind source separation may be viewed as a multiple-input, multiple- 
output (MIMO) blind channel estimation problem. 

Figure 3 is an illustration of a MIMO blind channel estimation scenario 
showing five unknown sources, Si, S2, S3, S4, S5, having distinct radiating patterns and 
five sensors, xi, X2, X3, X4, X5, having distinct receiving patterns. The sources, si, S2, S3, 

10 s 4 , S5, may provide and the sensors, Xi, X2, X3, X4, x 5 , may correspondingly receive, 
acoustic energy, electromagnetic energy, optic energy, mechanical energy, thermal 
energy, or a combination thereof. As shown in Figure 3, the five unknown sources, si, 
S2, S3, S4, S5, with distinct radiating patterns are generating a set of wavefields that are 
impinging on an array of five sensors, X|, X2, X3, X4, X5, with an unknown array 

15 manifold. Each source, si, S2, S3, S4, S5, provides a respective source signal. A BSS 
separation technique in accordance with the present invention, jointly extracts the set 
of source signals from an array of sensors (e.g., xi, X2, X3, X4, X5,) sampling the 
aggregate (composite) of the source signals' propagating wavefields at distinct spatial 
locations without knowledge of the signal characteristics or knowledge of the array's 

20 sensitivity as a function of direction of arrival or geometry. 

In order to develop a blind source separation technique suitable for separating 
narrowband signals given a set of outputs from an array of sensors with a relatively 
small spatial expanse and assess its performance, it is advantageous to develop a 
multiple-input multiple-output (MIMO) narrowband channel model for the array, state 
25 assumptions made, state the problem mathematically, and develop a set of measures to 
evaluate the technique. 

As such, a narrowband MIMO channel model is developed by starting with the 
most general convolutional MIMO channel model and then placing restrictions on the 
signal bandwidth and array size to simplify the problem, resulting in the narrowband 
30 model as utilized herein. Signal and noise assumptions are then presented and the 

blind source separation technique in accordance with the present invention is described 
mathematically and graphically. Two performance measures to be used in assessing 
the performance are then described including the novel concept of separation power 
efficiency (SPE). 
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Four multiple-input multiple-output (MIMO) channel models applicable to the 
blind source separation problem are described herein. These models are the general 
channel model, the non-dispersive direct path only channel model, the general finite 
impulse response (GFIR) channel model, and the narrowband channel model. The 
BSS technique in accordance with the present invention is then described utilizing the 
narrowband channel model. 

The General Channel Model : In the most general case, the output of each 
element is modeled as a summation of the M source signals each convolved with the 
impulse response of the channel between the output of the source and output of the 
sensor plus the additive Gaussian noise referenced to the sensors input. That is, 

M 

*»«-2>#('KW + *(') d) 

where * denotes convolution. The impulse response, v 0 (t) , of the channel between 

the output of the f h source and the i* sensor output may be time varying and account 

for such phenomena as multi-path propagation, dispersion, sensor time-varying 
response, source motion, sensor motion, etc. This can be written in matrix form as the 
general multiple input multiple output (MIMO) channel model 

x(/) = [x 1 (0 x 2 (t) ... x N (t)f 

v u(') - v im(0 



NM 



»* (0 



(2) 



= V(0*s(/)+n(0 



where [ ] denotes transposition. 



The Non-Dispersive. Direct Path Only Channel Model : When there is no 
multi-path, motion, or dispersion, the channel impulse response can be modeled by a 
delay and attenuation. That is, 



where a is the cascaded attenuation/gain from the output of f source to the i 



ith 



V 



sensor output and x i} is the propagation time (delay) from the output of j th source to 
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the output of the i th sensor. Under this model, when the sifting property of the delta 
function is employed, the output of the i th sensor (ignoring the noise) becomes 

*»(')=!># (')**, (') 

7=1 

-f,<h*('-*i) m 'j(t) W 

M 

At this point a "differential" delay is defined as the difference in propagation 
5 time from the output of the f h source to the output of the k th sensor and to the output 

of the I th sensor. 

AT u>j mT t-** ( 5 ) 

This differential time delay defines the time difference of arrival between two sensors 

for a given signal and is a measure of the spatial expanse of the array of sensors. 

10 Additionally, to facilitate situations when the minimum propagation delay from the 
/* source to the sensors is much greater than the maximum differential propagation 

delay, that is min(r (/ )l» max Ar lkj , the propagation time x i} is decomposed into 

two components, a "reference" delay, which is defined as the average propagation 
time from the output of the source to the output of the sensors and denoted as x } , and 

1 5 a "relative" delay, which is defined as the difference in propagation time between the 
reference time delay and the actual propagation time and denoted as Ax iy The 

propagation time from the /* source to the /'* sensor can then be expressed as 

T^Tj+ATy. (6) 

Figure 4 is a graphical illustration of time delays between sensors and sources. 

20 The decomposition of the propagation time as depicted in Figure 4 includes five 

sources, labeled s } , s 2 , . . . , s 5 , with associated reference delays fj , r 2 , . . . , r 5 , which are 

generating a set of wavefields that illuminate a set of five sensors, labeled x, , x 2 ,. . . , x s , 
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and the relative time delay, Ar^ , is shown for the first source, s x , and the third sensor, 
x 3 . Using the above definitions, the differential time delay can be reformulated as 

follows: 

= fc + 4r # )-(r y +4T # ) (7) 

5 Both the differential and relative time delays are utilized in the formulation of the 
narrowband and the general finite impulse response models. 

The General Finite Impulse Response (GFIR) Channel Model : The general 
model is often simplified by modeling the channel between the output of the f h 

source and the i' h sensor output, v fJ (t) , as a FIR filter or tapped delay line. As with 

10 the general model, the GFIR Model may be time varying and can account for such 

phenomena as multi-path propagation, dispersion, sensor time-varying response, 
system motion, etc. The FIR filter used to model v tJ (/) must be long enough to 

account for the multi-path delay spread of the channel as well as the relative time 
delay, Ar i , with a "reference" delay, r , accounted for by defining a delayed version 

15 of the source signal as it's input. That is the input to the set of FIR filters used to 
model the channels between the output of the j th source and array of sensors is 

0(')«*y('-*y)- ( g > 

The FIR filter or tapped delay line model is valid for a fading channel when the 
coherence bandwidth of such a channel is much less than the noise equivalent 

20 bandwidth of the source signal, that is BW NEq [s y (/)] < BW™ H , where the coherence 

bandwidth is defined as the reciprocal of the multi-path delay spread. In this situation 
the multi-path components in the channel separated by a delay of at least 

2nj BW NEq \jj (f )] are resolvable and the fading phenomenon is referred to as being 
"frequency selective". Thus the channel impulse response can be represented as 

25 v, (/) = £ vf (t)s(t - 2nljBW HEq [>, (/)]) (9) 

/=>0 

where the time varying complex- valued channel gain of the I th component can be 
represented as 

v«(0=«f(0^ (,) - do) 
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The length of the model, L y , is the number of resolvable multi-path components 
which is 

L v =\BW m [ Sj (t)]/BW? OH ] (11) 

where [ ] denotes the ceiling function. For the GFIR channel model, the length of 

5 the FIR filter has to not only accommodate the multi-path delay spread but also the 
relative time delay, Ar^ . That is equation (11) becomes 

h =[5^[5,(/)]{(|^|/2^) + (l/5^)]]. (12) 

In practice, the length of all the FIR filters are set to a common value, L , which is 
defined as 

L = max(L,.). (13) 

When the coherence bandwidth is greater than the noise equivalent bandwidth 
of the source signal, that is BW NEq [s y (/)] < BWf~ OH , the fading is referred to as 

"frequency non-selective" and the fading model reduces to a single time varying 
complex weight. That is L & =1 9 and thus 

v,(0=vj 0) (0 

= 4 0 )(/)^ 0,(,) (14) 



which begins to look like a time-varying narrowband model. However, for the above 
simplification to a single complex weight to hold in array signal processing, the source 
signal must have a noise equivalent bandwidth much less then the center frequency 
and the array of sensors must have a relatively small spatial expanse, that is 

20 B FF NEq [sj(t)] « ©j (15) 



max 



At 



^BW NEq [ Sj (t)\ (16) 



The Narrowband Channel Model : A measure of the spectral support of a 
signal is the noise equivalent bandwidth, denoted as BW NEq [ ] . By the duality 

principle of time and frequency, the inverse noise equivalent bandwidth can be used as 
25 a measure of the temporal support of the signal, in other words it is can be used as an 
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indication of the decorrelation time of the signal. When the signal noise equivalent 
bandwidth is much less then the center frequency, that is 

WNEq[Sj(t)] « COj (17) 

where co } is the center frequency of the /* source, then the propagation delay, or 

5 relative propagation delay, can be modeled as a phase shift. In this situation, when 
there is no dispersion or multi-path, the channel model is referred to as the narrowband 
model. 

However, since the phase shift is modulo 2n with respect to the center 
frequency, the requirement that the bandwidth be much less than the center frequency 
1 0 is itself insufficient in order for the time delay to be modeled as a phase shift and 
preserve the waveform, i.e. negligible inter-symbol interference (ISI) is induced in a 
digital communications signal. Therefore, for the narrowband model to hold, the array 
of sensors must also have a relatively small spatial expanse. That is 

AT likJ \« InlBW^ls^t)] (18) 



max 



1 5 Since Ar lkj = At 1} - At^ , requiring 



max \A Ty 



«*lBW NEq [ Sj (t)] (19) 



is a sufficient condition to guarantee (18) holds, via the triangle inequality. When the 
Narrowband conditions defined in (17) and (19) hold, the relative time delay is 
negligible in comparison to the decorrelation time of the signal and thus 

20 Sj ( t - tj - A t 9 ) = Sj ( f - tj ) (20) 

which says the waveform is preserved (within a phase shift). Thus, the relative time 
delay can be modeled as a phase shift, 

v<,( t )* s A t )= a A t - T >j)* 5 A t ) 

= a ljSj (t-T^ 

= a iJ s j (t-Tj-AT u ) (21) 

= aS*'j(t-fj)'W(t) 
where r j {t ) = Sj [t - tj \ , $ 9 = a>j A x i} , and a complex weight, v^. , is defined as 
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This complex weight together with the other N-\ weights associated with the j ,h 
signal form the /* steering vector. 

|T 



V 7=[ V U V 2j - ^j] 



(23) 



The output of the i* sensor is then 



M 



7=1 



(24) 



As done for the general m /* model, this can be re-formulated in matrix form for the 
vector of sensor outputs as 

x (0 = 1X0 M0 - MOT 
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(25) 



= Vr(0+n(0 

Due to conservation of energy, the total average signal power from the f h 
source illuminating the array can never exceed Pj . Since the signal-to-noise ratio is 

established at the input of the sensor, in the total array gain can be viewed as being 
normalized. Thus for the narrowband model, the inner product of the f h column of 

the mixing matrix V is, 

N 

H \p ♦ 



1=1 
N 



N 



(26) 



1=1 



= 1 



where [ J denotes the Hermitian transpose. 
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Signal and Noise Assumptions: The following assumptions are made about the 
source signals and noise vector. These assumptions are made to allow the use of the 
fourth-order cumulant and to ensure a sufficient number of degrees of freedom for the 
separation technique to exploit. Assumptions Al and A2 ensure that the fourth-order 

5 cumulant of the source signals exist. The zero-mean assumption is not necessary for 
the use of cumulants but is assumed since practical propagating electro-magnetic 
signals have a zero-mean. Assumptions A3 and A4 are particularly useful to the use 
of cumulants for the BSS problem. Without them, the noise sources would need to be 
treated as a signal source thus requiring additional degrees of freedom in the array. 

10 Note that the noise sources are not assumed to be temporally or spatially white. This 
is contrast to assumptions made in other second order techniques. The final 
assumption on the number of source signals helps to ensure that the there are enough 
degrees of freedom to perform the separation using the matrix-pencil approach. 

The first assumption (A\) : The M source signals illuminating the array are 
1 5 statistically independent non-Gaussian stationary random processes. Assumption Al 
is represented mathematically as follows. 

M 

f n ^..,r tl {r t ,r 2 . . .,r M ) = \{f rj (r, ) (27) 

The source signals are assumed stationary to order four over the estimation period of 
the spatial fourth-order cumulant matrix. 

20 The second assumption (AT) : The M source signals illuminating the array 

have zero mean with power Pj and a non-zero fourth-order moment. Assumption A2 



is represented mathematically as follows. 

£[m,(f)] = 0 (28) 

£[r, (0] = E^m, (t)] = ^[m, (/)] = 0 (29) 

25 £[m,.(r)m;(/)] = l (30) 

*[0 (00 (')] = E[tfmj UjPjrn) (/)] = J^[m y (t)m) (t)] = ? } (31) 

('K (<h (<K (<)] = P] E[ mj (t)m) (0 m, (t)m] (/)] (32) 

E[ mj {t) m] (t) mj (t)m) (f )] ^ 0 . (33) 
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The third assumption (A3) : The set of source signals (processes) and set of 
noise processes are statistically independent. Assumption A3 is represented 
mathematically as follows. 

M 

r„M n „('l>^> -->^"l>^ - ( 34 ) 

The fourth assumption (A4) : The noise processes are stationary zero-mean 
Gaussian random processes. They are not assumed to be spatially or temporally 
independent. Assumption A3 is mathematically represented as follows. 

n i (t)QN(Q,<rf) (35) 



■(0-h(0^(0v- >MO] t o*(o,k w ) 



fn,,n lr ..,n„ fa »**2v " * » W W ) — 



1 



(2;r) ,/2 |detK n 



|V2 



e -n T K„-'n 



(36) 



(37) 



The fifth assumption ( A5) : The number of sources is less than or equal to the 
number of sensors, i.e. M <N . 

Figure 5 is a functional block diagram of an apparatus for performing blind 

source separation (BSS) of M statistically independent narrowband source signals 

given a set of outputs from an array of sensors with a relatively small spatial expanse 

and with an arbitrary and unknown array geometry, in accordance with an embodiment 

of the present invention. The BSS technique as described herein determines a 

separation matrix W that will diagonalize the mixing matrix V . This involves 
finding a NxM separation matrix W , with complex elements w 0 , 



W = 



W 



NM J 



(38) 



that will diagonalize the mixing matrix, V . That is, a separation matrix W is desired 
such that the product W H V results in a MxM diagonal "loss" matrix with elements 
Pj- 
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W H V = 



Pi 
0 



0 

P2 



0 



0 



0 



0 

Pm 



(39) 



+ W H n(r) 



(40) 



When the separation matrix, W, is applied to the vector of sensor outputs, the 
result is 

y (t) = W H x(r) = W H {Vr + n(t)} 

= W H Vr(r) + W H n(r) 

p { 0 — 0 

0 p 2 \ : 

: \ \ 0 

0 - 0 p M \ 

+ W H n(/) 

and the source signals are separated. For mathematical completeness note that the 
vector r(t) e the vectors x(f),n(/) e C N , and the matrices V,We C NxM . If the 

loss matrix is the identity matrix, the separation process has captured all of the signal 
energy illuminating the array thus guaranteeing that the separated output signal has 
reached the maximum achievable signal to interference plus noise ratio. 

As developed previously in the narrowband model, the output of each sensor is 
a weighted linear mixture of the independent source signals plus noise. 



M 



*,(')=Xv/(') +"/(') 



(41) 



/=1 



Applying the separation matrix to the vector of sensor outputs separates the sources. 



= W H {Vr(r) + n(r)} 
= W H Vr(r) + W H n(r) 



(42) 



The / element of the separation process output vector, y, (r) , is an estimate of the 
/* source signal, r y (r) and is the inner product of the J* column of the separation 

matrix and the vector of sensor outputs. 
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(43) 



1=1 



Substituting equation (41) into equation (43) yields 



il/ [ M 

Z^lZviM + ^W 
I'-i 

M N 

= Z w </ Z v# (0 + Z ^ (0 



1=1 

AT M 



(44) 



i=l M 
N 



A> M 



f=l 



^ AT 

= Z v/j W + Z Z w l v * r i (') + Z w l n i (0 

r=l 



where it is clear there are three distinct terms corresponding to the desired signal, the 
residual interference, and the output noise. Of particular interest in evaluating the 
performance of communications and signal intelligence systems is the second-order 
moment of each of these terms. The second-order moment of the first term is the 
desired signal output power and is defined as 

2 

ZWy(0 



Sj = E 



N 



i=l 



= E 



r h 






• 






Zw)W 




li=l J 




U=i J 





(45) 



ZIwi¥[o (Oo (')] 

i=l A=l 



10 Applying assumptions Al and A2, 



(46) 



and thus equation (45) becomes 



N N 



(47) 



/=! A=l 



which can be represented using vector notation as 

The second-order moment of the second term in (44) is the residual 
interference power and is given by 



(48) 
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w 



J v? (0 



i=i /=i 



= E 



N M 




W (') 



N M 



v 



r=1 /=1 




1 m=l 
m*j 



(49) 



AT N 



= £ Z £ [ r / (0 r - W] Z Z w i w£» w * 



l*j m*j 



i=\ k=\ 



However, by assumption Al the signals are statistically independent and therefore 

ElnityZit)]^, formal. (50) 

Additionally, applying the stationarity of assumption Al and assumption A2, 

E[ r >(t)r;(t)]^Pr (51) 

Using (50) and substituting equation (51) into equation (49), the residual interference 
power reduces to 

M N N 
/=! /=1 k=\ 



which can be represented using vector notation as 

M 



(53) 



/=1 



The second-order moment of the third term in (44) is the output noise power 
and is given by 



Nj=E 



Z^(') 



N 



N 



(54) 



N N 



i=i *=i 
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which can be represented using vector notation as 

^=w 7 H £[n(/)n H (/)]w y .. (55) 

By definition and assumption A4, the expectation of the outer product of the noise 
vector is the noise covariance matrix, 

5 E[n(t)n" (t)] = K„ (56) 

and thus the output noise power is 

A^wjK^w.. (57) 

To evaluate the effectiveness of a blind source separation technique, a measure 
of the quality of the separation is utilized. As previously described, the blind source 
1 0 separation technique as described herein determines a separation matrix W that will 
diagonalize the mixing matrix V . Two measures to assess the quality of a blind 
source separation algorithm are developed herein. Performance of the BSS technique 
may be measured in terms of residual interference and in terms of the efficiency of the 
algorithm in "capturing" all available signal power illuminating the array of sensors. 

1 5 One measure of the quality of separation is the amount of residual interference 

found in a signal output after the separation matrix has been applied. Specifically, the 
power of the residual interference relative to the desired signal in the estimate of the 
desired source averaged over all sources as well as the peak or maximum residual 
interference-to-signal ratio to assess the separation technique in terms of suppressing 

20 co-channel interference are proposed for use. This measure is of significance because, 
if the separation matrix does not perfectly diagonalize the mixing matrix, the off 
diagonal terms of the resultant matrix will permit residual interference in the signal 
outputs. 

In most communications applications, the common measure of the amount of 
25 interference is the signal-to-interference ratio, which is the ratio of the desired signal 
power to the combined power of all interfering signals. However, as the goal of the 
blind source separation is to completely eliminate all interference, this ratio could 
become extremely large. As a result, the Interference-to-Signal ratio (ISR), which 
quantifies the residual power of the interference that a blind source separation 
30 algorithm or technique fails to suppress relative to a particular desired signal power, is 

* 

proposed. The better an algorithm is the smaller this ratio will become. 
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The ISR of a particular desired signal is defined as 

< 58 > 

j 

Substituting (53) and (48) into (58), the ISR for a particular signal is 

M 

Sj = '\ H h • (59) 

Pw - V V w . 

J J J J J 

5 This value is also known as the rejection rate. 

The overall quality of separation of a blind source separation technique may be 
measured by looking at the mean value of the individual source signal ISR's, £ } , over 

all j . Thus the primary measure to be used to evaluate the performance of a blind 
source separation algorithm in terms of residual interference will be the average ISR 
1 0 given by 

The secondary measure in terms of residual interference will be the peak or 
maximum ISR, which is defined as 

ISR m „=max[Cj]. (61) 

15 This secondary measure ensures that all source signals are effectively separated with 
an ISR no worst than ISR^ . 

A second measure of the quality of separation is utilized to determine the 
efficiency of the source separation matrix in terms of its ability to make use of the 
available signal power. A BSS technique is considered more efficient if the output 
20 signal-to- interference-plus-noise ratio is maximized, thus having greater sensitivity in 
terms of being able to capture smaller signals, than a BSS technique not maximizing 
the output signal-to- interference-plus-noise ratio. 

The efficiency of a blind source separation algorithm in using all of a source's 
signal power illuminating the array of sensors is yet another important measure of its 
25 quality of separation. This measure determines how much of the available signal 
power from a particular source is wasted or lost in the separation process. This loss 
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results in a lower signal-to-noise-plus-interference ratio then would otherwise be 
theoretically achievable and thus a loss in system sensitivity. The Separation Power 
Efficiency (SPE) for a particular separation process output relative to the desired 
source signal's available normalized power is defined as 

5 £ yS A where (62) 

f j 

tj is indicative of the separation power efficiency for the j 01 source of the plurality of 
sources, Sj is indicative of a power of a separated signal from the f 1 source, and Pj is 
indicative of a normalize power of a signal from the 3 th source. 

Substituting equation (48) in for the separation process output power reveals that the 
10 particular SPE 

r _ J J J J TT J 

Pj (63) 

II ii 

= w V V w 

J J J J 

depends only on the steering vector for the j th source and the j th column of the 

separation matrix. As with ISR, both the average SPE and the minimum SPE, defined 
as 

is Vifl'' (64) 



and 



^min sm J n fe] (65) 



respectively, will be used to evaluate the separation power efficiency. 

Note that by the definition of the illuminating source signal power, i> , that the 

20 maximum value the SPE can achieve is one. Thus the maximum achievable average 
SPE is also one. A separation algorithm that achieves an SPE of one is guaranteed to 
have maximized the source signal power in the corresponding separation process 
output. The minimum SPE provides a measure of ensuring that all sources are 
successfully separated with a minimum separation power efficiency. 

25 A BSS technique in accordance with an embodiment of the present invention 

utilizes cumulants, specifically spatial fourth order cumulant matrices. To better 
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20 



25 



understand the use of cumulants in performing blind source separation, a cumulant 
definition and associated properties are provided below. 

The joint cumulant, also known as a semi-invariant, of order N of the set of 
random variables {s } , s 2 , . . . , s N } is defined as the N th -order coefficient of the Taylor 

series expansion about the origin of the second characteristic function. See, for 
example, C.L. Nikias and A.P. Petropulu, Higher-Order Spectra Analysis: A Non- 
Linear Signal Processing Framework. (PTR Prentice-Hall, Upper Saddle River, NJ: 
1993) and M. Rosenblatt, Stationary Sequences and Random Fields (Birkhauser, 
Boston, MA: 1985), which are hereby incorporated by reference in their entirety as if 
presented herein. The second characteristic function is defined as the natural 
logarithm of the characteristic function, 

^,(tf>i>*>2>- ..,^) = ln[0 8 (6? 1 ,6; 2 ,.. . 7 a> N )] (66) 



where the characteristic function is defined as 



(67) 



1 5 The joint N' h -order cumulant is then 



Cum[s l ,s 2 ,...,s N \ = {-j) - 



dco x do) 2 --d(o N 



(68) 



Cumulants, unlike moments, cannot be directly estimated from the data. See, 
for example, A.K. Nandi, Blind Estimation Using Higher-Order Statistics (Kluwer 
Academic, Dordecht, The Netherlands: 1999), which is hereby incorporated by 
reference in its entirety as if presented herein. However, cumulants can be found 
through their relationship to moments and can thus be estimated indirectly by first 
estimating the required moments. The relationship of cumulants to moments is 
described in M. Rosenblatt, Stationary Sequences and Random Fields (Birkhauser, 
Boston, MA: 1985) for the N' h -order joint cumulant of the set of random variables 



N 



Hp) p 



c«m[,„5 2 ^]=K-ir(p-0'iin £ 



(69) 
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where there are N(p) ways of partitioning the set of integers 
{l, 2,. . . , N) into p groups, each denoted as g t p n , such that 

fW={ 0 } 

';' (70) 

As an example, for the case N = 4 , the partitioning is defined on the set of 
5 integers {1,2,3,4} and is given in Table 1.0 below. 

Table 1.0: All Possible Partitions for N = 4 



p 


Hp) 


oM:p t p/ial;N[p) 


1 


i 


{1,2,3,4} 






{1}{2,3,4};{2}{1,3,4};{3}{1,2,4};{4}{1,2,3}; 


2 


7 


{1,2}{3,4};{1,3}{2,4};{1,4}{2,3} 






{1}{2}{3,4};{1}{3}{2,4};{1}{4}{2,3}; 


3 


6 


{2}{3}{1,4};{2}{4}{1,3};{3}{4}{1,2} 


4 


1 


W{2}{3}{4} 



The 4 -order joint cumulant as a function of the moments is then 

Cum [ j, , s 2 , s 3 , s A ] = E [5,j 2 s 3 s A ] - E [s, ] • E fas^ ] - E [ s 2 ] • E [s^s 3 s 4 ] 

-E [s 3 ]E [s i s 2 s 4 ] - E[s A ] ■ ] - E [s^ ] • E [ j 3 $ 4 ] 

-E [ v 3 ] E [s 2 s 4 ] - E [s t s 4 ] ■ E [s 2 s 3 ] + 2E [s t s 2 ] ■ E [s 3 ] ■ E [s 4 ] 

+2E [ v 3 ] • E [s 2 ] ■ E [s 4 ] + 2E [ Sl s 4 ] ■ E [s 2 ] • E [s 3 ] 

+2E [s 2 s 3 ] • E [s, ] ■ E [s 4 ] + 2E [s 2 s 4 ] • E ] • E [s 3 ] 

+2E[s 3 s 4 ]E[s l ]E[s 2 ]-6E[s i ]E[s 2 ]E[s 3 ]E[s 4 ] 

1 0 Note that equation (71) shows that computation of the N -order joint cumulant 
requires knowledge of all moments up to order N . 

Cumulants possess several properties that make them attractive for use in the 
blind separation of a linear mixture of unknown statistically independent signals in 
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spatially and/or temporally correlated Gaussian noise, especially at a low signal-to- 
noise ratio. 

One property that makes cumulants attractive for use in blind source separation 
is that if the set of random variables { s v s 2 , . . . , s N } can be divided in to two or more 

5 groups that are statistically independent, then their N ih -order joint cumulant is zero. 
Thus, the cumulant operator in the blind separation of statistically independent sources 
will suppress all cross-source signal cumulant terms. In general, this is not the case for 
higher-order moments. Another property that makes cumulants attractive for use 
inBSS is that the Cum[s x + n x ,s 2 + h 2 ,...,s a , = 

1 0 Cum [s { , s 2 , . . . , s N ] + Cum [n, , ^ , . . . , n N ] . Because in general the set of signal terms 

{s x , s 2 , . . . , s N } and the set of noise terms {n, , ^ , . . . , n N } are statistically independent 

from each other, the N th -order joint cumulant of the terms of their vector sum, 

{s } + n x , s 2 + , . . . , s N + n N } , is the sum of their individual joint cumulants. Therefore, 

the cross cumulants between the noise terms and signal terms will be zero. This 
1 5 property is important in guaranteeing that the spatial fourth-order cumulant matrix can 
be decomposed into the sum of two matrices, one corresponding to the signals and the 
other corresponding to noise vector. 

Yet another property that makes cumulants attractive for use in BSS is that the 
joint cumulant of order N>2 of a Gaussian random variable is zero. Because the noise 
20 vector is a multi-variate Gaussian random process, hi = f n x 

its joint cumulant of order three or higher will be zero. That is Cum [«, , , . . . , n N ] = 0 . 

This last property results in the spatial fourth-order cumulant matrix not having a noise 
subspace and the only non-zero elements of the matrix are associated with and only 
with the source signals. This is true even if the noise vector is spatially or temporally 
25 correlated. 

Finally, cumulants of order higher than two preserve phase information that is 
lost by the use of second-order statistics, such as correlation. For example, auto- 
correlation destroys the information necessary to distinguish between minimum phase 
and non-minimum phase signals. Thus, two signals may have identical second-order 
30 statistics yet have different higher-order statistics. This property is of particular 
interest in handling signals with identical auto-correlation functions and adds 
additional degrees of freedom for finding a set of time lags where a group of source 
signals will have different higher-order cumulants. This property is particularly 
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advantageous to a BSS technique in accordance with the present invention because a 
condition of identifiability of this BSS technique is that all signals have a unique 
normalized fourth-order auto-cumulant. Note that the fourth-order cumulant is used 
because odd order cumulants of a process with a symmetric distribution will be zero. 

5 Four properties of cumulants utilized in the BSS technique in accordance with 

the present invention are described below. Proofs of these cumulant properties may be 
found in C.L. Nikias and A.P. Petropulu, Higher-Order Spectra Analysis: A Non- 
Linear Signal Processing Framework. (PTR Prentice-Hall, Upper Saddle River, NJ: 
1993) and M. Rosenblatt, Stationary Sequences and Random Fields (Birkhauser, 

10 Boston, MA: 1985). 

Cumulant Property 1 : 

The N th order joint cumulant of the set of random variables [a x s x , a 2 s 2 , . . . , a N s N } is 

f N 1 

Cum[a { s x ,a 2 s 2 ,... :> a N s N ] = \Y[ a i \Cum[s }9 s 29 ... 9 s N ] where {a l9 a 29 ... 9 a N } are 



constants. 



15 Cumulant Property 2 : 

If the set of random variables {j„j 2> — iff*} can be divided in to two or more 

groups that are statistically independent, then their N th -order joint cumulant is 
zero. 

Cumulant Property 3 : 

20 If the sets of random variables {s X9 s 29 ... 9 s N } and {n,,^,...,^} are statistically 
independent, i.e. f„(s l9 s 29 ... 9 s„ 9 n i9 n l9 ... 9 n N ) = f B (s l9 s 2 ^ 
then the AT" 1 -order joint cumulant of the pair-wise sum is 
Cum[s } + n^ 9 s 2 + n 29 . . + n N ] = Cum [ s t , s 2 , . . . , s N ] + Cum [ n, , n 2 , . . . , n N ] . 

Cumulant Property 4 : 

25 If the set of random variables {w,,^ n N ) are jointly Gaussian, then the joint 

cumulants of order N > 2 are identically zero. That is, if 

n = [n l ,n 2 ,...,n N ] T U N(|i„,K n ) , then Cum[n l ,n 2 ,...,n N ] = 0. 
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A BSS technique in accordance with the present invention utilizes a fourth 
order spatial cumulant matrix. Three definitions of the spatial fourth-order cumulant 
matrix and associated properties are provided below. 

The spatial fourth-order cumulant matrix is used as a basis for estimating a 
5 separation matrix at low signal-to-noise ratios and in the presence of spatially and 
temporally correlated noise since it theoretically has no noise subspace, even if the 
noise is correlated. This eliminates the need to use either degrees of freedom and/or 
secondary sensor data to estimate the noise subspace, which must be removed in order 
for the matrix-pencil to be formed. As described below, the absence of the noise 
10 subspace is a direct result of using a higher-order cumulant, i.e. order > 2, and is 
particularly advantageous to a blind source separation technique in accordance with 
the present invention. 

The three spatial fourth-order cumulant matrix definitions and their properties 
are presented herein with consideration of the fact that the sensors are in reality never 

1 5 omni-directional, never have identical manifolds, and that different sets of time lags 
are needed to estimate a pair of spatial fourth-order cumulant matrices to form the 
matrix-pencil. These considerations are a clear distinction from previous treatments of 
the spatial fourth-order cumulant matrix. See, for example, H.H. Chiang and C.L. 
Nikias, "The ESPRIT Algorithm with Higher-Order Statistics," Proc. Workshop on 

20 Higher-Order Spectral Analysis, Vail, CO., Jun. 1989, pp. 163-168, C.L. Nikias, C.L. 
Nikias and A.P. Petropulu, Higher-Order Spectra Analysis: A Non-Linear Signal 
Processing Framework (PTR Prentice-Hall, Upper Saddle River, NJ: 1993), M.C. 
Dogan and J.M. Mendel, "Applications of Cumulants to Array Processing - Part I: 
Aperture Extension and Array Calibration," IEEE Trans, Signal Processing, Vol. 43, 

25 No. 5, May 1995, pp. 1200-1216, and N. Yuen and B. Friedlander, "Asymptotic 
Performance Analysis of ESPRIT, Higher-order ESPRIT, and Virtual ESPRIT 
Algorithms," IEEE Trans. Signal Processing, Vol. 44, No. 10, Oct. 1996, pp. 2537- 
2550. Understanding the properties of the spatial fourth-order cumulant matrix such 
as its rank, null spaces, etc., and its relationship to the mixing matrix are beneficial to 

30 developing a signal subspace blind separation technique using fourth-order cumulants 
and a matrix-pencil in accordance with the present invention. 

A brief review of the spatial correlation matrix and its properties are provided 
below to aid in understand its use in a BSS technique in accordance with the present 
invention. The spatial correlation matrix of the sensor array output is defined in D.H. 
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Johnson and D.E. Dudgeon, Array Signal Processing: Concepts and Techniques. 
(PTR Prentice-Hall, Englewood Cliffs, NJ: 1993), which is hereby incorporated by 
reference in its entirety as if presented herein, as: 

R x (T) = E[x(t)x»(t-T)] (72) 

5 Substituting (25) for x(/) in to equation (72) and applying assumptions Al and A3, the 
spatial correlation matrix becomes 

R I (r) = £[{Vr(r) + n(r)}{Vr(r-r) + n(/-r)} H ] 

= £[Vr(/)r H (/-r)V H ]+ E[\r(t)n H (t-r)] 
+E[n(t)r H (t-T)\ H ] + E[n{t)n H (t-T)] (73) 
= VE [r (t ) r H (/ - t )] V H + E [n (t) n H (t - r)] 
= VR r (r)V H +R n (r) 

which has elements 

fr ML = £ v * v A r J (') 0 (' " 0] + E (') < (f - r)] (74) 

10 where the subscript rc indicates the element is in the r th row and c lh column. Since 
the signal and noise processes are assumed to be zero mean, assumptions A2 and A4, 
the spatial correlation matrix defined in equation (72) is equivalent to the spatial 
covariance matrix, and thus the terms are used interchangeably. 

In general, most second-order techniques make use of the spatial correlation or 
1 5 covariance matrix only at a delay lag of zero, {r = 0} . In such a case the spatial 

correlation matrix is Hermitian and non-negative definite. See for example D.H. 
Johnson and D.E. Dudgeon, Array Signal Processing: Concepts and Techniques. 
(PTR Prentice-Hall, Englewood Cliffs, NJ: 1993), C.L. Nikias and A.P. Petropulu, 
Higher-Order Spectra Analysis: A Non-Linear Signal Processing Framework. (PTR 
20 Prentice-Hall, Upper Saddle River, NJ: 1993), and A. Papoulis, Probability, Random 
Variables, and Stochastic Processes. (WCB/McGraw-Hill, Boston, MA: 1991), for 
example. Further, if the sensor outputs are linearly independent, that is 
E [{a T x (/)} {a T x (t )} * J > 0 for any a = [a, , a 2 , . . . , a N f * 0 , then the spatial correlation 

matrix is positive definite. As a consequence of the spatial correlation matrix being 
25 non-negative definite for r = 0 , its determinant will be real and non-negative, and will 
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10 



15 



be strictly positive if and only if the sensor outputs are linearly independent. 
However, if r * 0 then the spatial covariance matrix is indefinite and non-Hermitian. 

Spatial Fourth-Order Cumulant Matrix Definition 1 

The first definition of a spatial fourth-order cumulant matrix presented takes advantage 

of the steering vectors having a norm of one. This is stated mathematically in equation 

(26). As will be shown, this is utilized to factor the spatial fourth-order cumulant 

matrix into Hermitian form when the sensors are not omni-directional with identical 

manifolds. The first spatial fourth-order cumulant matrix is defined at the set of time 
lags (r,,r 2 ,r 3 )as 



(75) 



i=l 



and is referred to as spatial fourth-order cumulant matrix 1. 



The spatial fourth-order cumulant matrix 1 as defined in (75) is in general a 
complex NxN matrix with the element in the r* row and c th column given by 



[C 4 x (7- | ,T 2 ,i-3)] re =X Cw7M [^( / - r i) X '('- T 2)^(0 JC c('- T 3)j 

where { }* denotes complex conjugation. Substituting equation (24) into (76), 
element rc becomes 



(76) 



N 

[c^(r„r 2 ,r 3 )] re = X C "' M 



ZVj (')+*(') 



r m y 

\m=\ ) 



f M V 

XvU*-*i) +, %( / - r i) 



v/=i 



(77) 



Then, by Cumulant Property 3 and assumption A3, (77) becomes 



[c:(^ 2 ,r 3 )] re =£C M m 



M 



M 



M 



M 



S V r/j WZ (' T l )Z V // r / (' " T 2 ) 

7=1 



Z V ^(^ T 3) 



ro=l 



(78) 



+ £ Cum (f) < (/ - r, ) n ( . (/ - r 2 r, )] 



37 



kJ 



WO 03/090127 



PCT/US03/11473 



where the terms have been re-ordered. However, by assumption A4 and Cumulant 
Property 4, 

2C M m[/i r (0«;(/-r 1 )« i (^-r 2 )<(/-r 3 )] = 0 (79) 

1=1 

and thus (78) becomes 

[C 4 , (r, , r 2 , r 3 )] = £ Cum \ £ V , (f )£ vtf (/ - r, )£ Vl (f - 1 2 ) 

M U» W M (g0) 



AY 1 

!>;„,/•; (/-r 3 ) 



Then, by the source signals statistical independence of assumption Al and repeatedly 
applying Cumulant Property 3, equation (80) reduces to 



Using Cumulant Property 1 , the complex weights may then be pulled out in 
1 0 front of the cumulant operator in equation (8 1 ) to give 

[c:(r,,r 2 ,r 3 )] re =zf>^ (82) 



Reordering the summation yields 

M N 

[cl(r u T 2 ,T 3 )] n =Zv,v;IvXC t /m[r ; (/)r;(/-r,)r y (/-r 2 )r;(r-r 3 )]. (83) 

7=1 '=1 

N N 

However, since the steering vectors have a norm of 1, that is ^v*v l7 = J/r£ - 1 » 

/"=! /=! 

1 5 equation (83) reduces 

[c:(r M r 2 ,r 3 )] =Zv,v;.C W m[r / (/)r;(/-r 1 )r / (/-r 2 )r,-(/-r 3 )]. (84) 

From (84) it can be seen that spatial fourth-order cumulant matrix 1 can be 
factored into Hermitian form, as was the case for spatial correlation matrix, 

C:(r„r 2 ,r 3 ) = VC:(r 1 ,r 2 ,r 3 )V H (85) 

20 
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where C*(r,,r 2 ,r 3 ) is a MxM diagonal matrix with elements, 

[c; (t iM )] j = Cum [ rj (t) r] (t - x x ) r } (t - r 2 ) r/ (f - r 3 )] 

scJ(T„r 2 ,r 3 ) 

Expanding equation (85) it is found that spatial fourth-order cumulant matrix 1 can be 
written as a sum of the steering vector outer products scaled by the individual source 
5 signal's fourth-order cumulant. 

M 

C 4 X (*i > *2 . *i ) = Z c n fa > T * > r 3 ) v j v " 5 where (87) 

C\ (r t , r 2 , r 3 ) is the spatial fourth order cumulant matrix having a first time lag, xi , a 

second time lag, x 2> and a third time lag, x 3 , each time lag being indicative of a time 

delay from one of the plurality of sources to one of the plurality of elements; M is 
10 indicative of a number of sources in the plurality of sources; c* (r,, r 2 ,r 3 ) is a fourth 

order cumulant of a j th source signal from one of the plurality of sources having delay 
lags xi, X2, and X3; and v^v" is indicative of an outer product of a f 1 steering vector. 



From equation (87) it is clear that spatial fourth-order cumulant matrix 1 lies in the 
signal subspace spanned by the set of steering vectors. Note that the spatial fourth- 
15 order cumulant matrix does not have the noise subspace that is present in the spatial 
correlation matrix. What was the noise subspace in the spatial covariance matrix is 
now the nullspace of the spatial fourth-order cumulant matrix. This property will be 
shown to be true for the other spatial fourth-order cumulant matrix definitions 
presented. 

20 Spatial Fourth-Order Cumulant Matrix 1 Properties 

Spatial fourth-order cumulant matrix 1, C* (r,,r 2 ,r 3 ) , has several properties, an 

understanding of which will facilitate the development of a method for estimating a 
separation matrix W . Establishing the spatial fourth-order cumulant matrix l's 
matrix properties is a first step to the use of the generalized eigen decomposition of the 
25 matrix-pencil formed by a pair of spatial fourth-order cumulant matrix l's at two sets 
of time lags. Such things as its rank and its subspaces relationships to the mixing 
matrix's subspaces are advantageous in developing a signal subspace separation 
algorithm. Particular attention is paid to the fact the individual sensors are not 
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assumed to be omni-directional with identical directivity for each impinging source 
signal wavefield. 

Property 1 : Spatial fourth-order cumulant matrix 1 is Hennitian if and only if 
r x =r 2 -r andr 3 =0, i.e. C*(r,r,0). 

5 Property 2 : The trace of spatial fourth-order cumulant matrix 1 equals the sum 

of the signal fourth-order cumulants, which is the trace of the diagonal matrix 

M 

V j=i (88) 

= /r(c:(r 1 ,r 2 ,r 3 )) 

Property 3 : The column space of spatial fourth-order cumulant matrix 1 , 
1 0 denoted as C (C 4 X (r x , t 2 , r 3 )) , is spanned by the set of steering vectors. 

^(Cfc:^,^,^)^^^,,...^,,} (89) 

Further, if the mixing matrix has full column rank, then the set of steering vectors are 
linearly independent and they form a basis for the column space of spatial fourth-order 
cumulant matrix 1 . 

1 5 Property 4 : If V has full column rank, then the rank of spatial fourth-order 

cumulant matrix 1 equals the rank of the mixing matrix. That is 

p(c:(r„r 2 ,T 3 )) = p(\) (90) 

ifp(V) = Af , where p( ) denotes rank. 

Property 5 : The "right" nullspace of spatial fourth-order cumulant matrix 1 
20 and the "left" nullspace of the mixing matrix are equal if the mixing matrix has full 
column rank. 

N ^(r^r^N ,(V) (91) 

Spatial Fourth-Order Cumulant Matrix Definition 2 

The second definition for a spatial fourth fourth-order cumulant matrix is one 
25 modified from the definition described in H.H. Chiang and C.L. Nikias, "The ESPRIT 
Algorithm with Higher-Order Statistics," Proc. Workshop on Higher-Order Spectral 
Analysis, Vail, CO., Jun. 1989, pp. 163-168 and C.L. Nikias and A.P. Petropulu, 
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Higher-Order Spectra Analysis: A Non-Linear Signal Processing Framework. (PTR 

Prentice-Hall, Upper Saddle River, NJ: 1993). These definitions are used and the set 
of time lags (r,,r 2 ,r 3 ) are incorporated to obtain spatial fourth-order cumulant 

matrix 2. 

cf(r 19 r 2) r 3 ) S C^[{ X (f)x (t-T x )x(t-T 2 )}x H (t~T,)] (92) 

Spatial fourth-order cumulant matrix 2 is slNxN matrix with the element in 
the r' h row and c th column 

[tf (*.,*,,*,)] =Q«i[x r (/)x;(l-r 1 )x r (f-r 2 )x;(f-r,)] . (93) 
Substituting equation (24) for x, (t) in equation (93), element rc becomes 



f 



M 



M 



V 



Ev>( r )+^(0 lE v * r *('- r 0 + ^('-*i) 

KM 



M 

/=1 J\ m=\ 



V M 



v 



j 



(94) 



Following the simplification of spatial fourth-order cumulant matrix 1 , Cumulant 
Property 3 and assumption A3 are applied to reduce equation (94). 



Cj'(r lf r 2 ,r 3 ) 



-Cum 



-ire 



M 



M 



M 



E v, WE v >* (< - *\ )E v, (' - r > ) 



A/ 



E v « r -('"" r 3) 



(95) 



+Ciiw[« r (/)«;(f-T 1 )» r (/-r 2 )n;(r-T 3 )] 



However, by assumption A4 and Cumulant Property 4, 
15 Cum[n r {t)n;(t-r,)n r {t-r 2 )n:(t- h )] = 0 

and thus (95) reduces to 



(96) 



c i'(*i»*2»*3) =Cum 

Jrc 



M 



ivi(of;«(/-fi)Zvi('-»a) 

_ >1 



z 



M 



(97) 
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Then, by the statistical independence of the source signals of assumption Al and 
repeatedly applying Cumulant Property 3, equation (97) reduces to 

C 4 ,' (r„r 2 ,r 3 )l = ffimU^ (t^rj (f-r, ) V , {t-r 2 )v c/ ; (t-rj]. (98) 

L Arc >1 L 

Using Cumulant Property 1 , the complex weights may then be pulled out in 
front of the cumulant operator in equation (98) to give 

fe' (r p r 2 , r 3 )1 = £ v.v^Cum [r, (/) r] (/ - r, ) r ) {t - r 2 ) r] (t - r 3 )] . 

L Arc j=\ 

However, v v = a and equation (99) reduces 



(99) 



|c 4 ;(r„r 2 ,r 3 )] =f i alv d v cj Cum[r J (t)rJ (/- Tl )r.(r-r 2 )r;(/-r 3 )] 

L Arc j=l 



(100) 



From (100) it can be seen that spatial fourth-order cumulant matrix 2 in 
general can not be factored into Hermitian form, as was the case for spatial fourth- 
order cumulant matrix 1 and the spatial covariance matrix. However, if 



v„. = a H v ri 

KJ KJ KJ 



(101) 



is defined, it can be factored in to bilinear form. 

C:'(^,r 2 ,r3) = VC:(r I ,r 2) r 3 )V H 



(102) 



where the element in the r' h row and c' h column of the NxM "modified" mixing 
matrix V is 



(103) 



Expanding equation (102), it is found that spatial fourth-order cumulant matrix 2 can 
be written as a sum of the outer products of the "modified" steering vector, v. , and 

steering vector scaled by the individual source signal's fourth-order cumulant. 

c: f (r 1 ,r 2 ,r 3 ) = X<(^^2^3)v 7 v7 0° 4 ) 
Note that the "modified" steering vector v^. is the j lh column of the matrix V . 

A question pertaining to spatial fourth-order cumulant matrix 2 is whether or 
not it is rank deficient. Following the derivation of the rank of spatial fourth-order 
cumulant matrix 1 , the rank of spatial fourth-order cumulant matrix 2 will be equal to 
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the rank of the mixing matrix if "modified" mixing matrix, V , and the mixing matrix 
both have full column rank. The mixing matrix V can be assumed to have full 
column rank since this can be guaranteed by design of the array. However, the rank 
of V cannot be guaranteed by design and as of yet, it is unclear if guaranteeing that 
5 the mixing matrix has full column rank is sufficient to guarantee that the "modified" 
mixing matrix will have full column rank. Although the "modified" mixing 
matrix V is the Hadamard product 

Vs V© V© V (105) 

the rank of the mixing matrix is not necessarily preserved. See for example, J. R. 

10 Schott, Matrix Analysis for Statistics. (John Wiley and Sons, New York, NY: 1997). 
At this point it shall be assumed that the Hadamard product preserves the rank of the 
mixing matrix and therefore that the mixing matrix having full column rank is 
sufficient to guarantee that the "modified" mixing matrix has full column rank. The 
implications of the "modified" mixing matrix not having full column rank will be 

1 5 clear in the subsequent sections. 

If the "modified" mixing matrix has full column rank, by inspection of 
equation (104) it is obvious that spatial fourth-order cumulant matrix 2 lies in the 
signal subspace spanned by the set of "modified" steering vectors. Again, the noise 
subspace in the spatial covariance matrix is now a nullspace of spatial fourth-order 

20 cumulant matrix 2. Note that in H.H. Chiang and C.L. Nikias, "The ESPRIT 

Algorithm with Higher-Order Statistics " Proc. Workshop on Higher-Order Spectral 
Analysis, Vail, CO., Jun. 1989, pp. 163-168 and C.L. Nikias and A.P. Petropulu, 
Higher-Order Spectra Analysis: A Non-Linear Signal Processing Framework. (PTR 
Prentice-Hall, Upper Saddle River, NJ: 1993), the elements/sensors are omni- 

25 directional with unit gain so that a] = 1 and as such, spatial fourth-order cumulant 

matrix 2 and spatial fourth-order cumulant matrix 1 would be equal and the 
"modified" mixing matrix has full column rank. However, this is an unrealistic 
assumption since in practice sensors are never omni-directionai. 

Spatial Fourth-Order Cumulant Matrix 2 Properties 

30 If the "modified" mixing matrix V has full column rank, spatial fourth-order 

cumulant matrix 2 will possess many of the same properties that spatial fourth-order 
cumulant matrix 1 does. The subsequent sections derive the key properties associated 
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with the development of a matrix-pencil signal subspace separation technique with the 
assumption that the "modified" mixing matrix has full column rank. 

Property 1 

Spatial fourth-order cumulant matrix 2 is in general non-Hermitian. It will be 
5 Hermitian if and only if r, = r 2 = r andr 3 = 0 , i.e. C* (r,r,0) and the sensors all 
have the identical gain for a given signal. 

Property 2 

The trace of spatial fourth-order cumulant matrix 2 equals the sum of the signal 
fourth-order cumulants scaled by the sum of the sensor magnitudes to the fourth 
10 power. 



tr 



f \ M N 

< (r„r 2 ,r 3 ) (r,,r 2 ,r 3 ) (106) 

V / j=] r=\ 



Property 3 

The column space of spatial fourth-order cumulant matrix 2, denoted as 
C ^cf (r,,r 2 ,r 3 )j , is spanned by the set of "modified" steering vectors. 



15 sp 



C^(r p r 2 ,r3)J = {v l5 v 2 ,...,v M } (107) 



Further, if the "modified" mixing matrix has full column rank, then the set of 
"modified" steering vectors are linearly independent and they form a basis for the 
column space of spatial fourth-order cumulant matrix 2. 

Property 4 

20 The rank of spatial fourth-order cumulant matrix 2 equals the rank of the mixing 
matrix, if V and V have full column rank. That is 

p(cf(r p r 2 ,r 3 )j = p(V) (108) 
if p(V) = p(v) = M , where p( ) denotes rank. 
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Property 5 

The "right" nullspace of spatial fourth-order cumulant matrix 2 and the "left" 
nullspace of the mixing matrix are equal if the mixing matrix and "modified" mixing 
matrix have full column rank. 



N 



C 4 ;(r l5 r 2 ,r 3 )j = AT / (V) 



(109) 



Spatial Fourth-Order Cumulant Matrix Definition 3 

The third and final definition for a spatial fourth-order cumulant matrix 
incorporates the time lags (r p r 2 ,r 3 ) and results in the following equation. 

cf (r 1 ,r 2 ,r 3 )^C W /n[{x(r)x (t-r x )x (,-r 2 )}x T (f-r 3 )] 



(110) 



Spatial fourth-order cumulant matrix 3 is again a NxN matrix with the element 
in the r th row and c lh column 



[cf (r p r 2 ,r 3 )l =Cww[x r (/)x;(r^r I )x;(r-r 2 )x c (/-r 3 )]. 

L Arc 

Substituting equation (81) for *,(*) in equation (111), element rc becomes 
c f { r \^h) = Cum 

L Arc 



(111) 



f M 

Z V ^/( r - T 2) + ^(^^) 



V / M 



\ 



(112) 



Following the simplification of spatial fourth-order cumulant matrix 2, Cumulant 
Property 3 and assumption A3 are applied to reduce equation (1 12). 



Cf (t x ,t 29 Ti) 



-Cum 



-J re 



M 



M 



M 



/=1 



A/ 



iw=l 



(113) 



+Cum[n r (t)n r {t-r x )n r (t-r 2 ) n c (t - r 3 )] 

However, by assumption A4 and Cumulant Property 4, 

C«m[/i f (r)/i;(/-r 1 )»;(r-r 2 )« c (r-r 3 )] = 0 (1.1) 
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and thus (95) reduces to 



[cf (*i,r 2 ,r 3 )J =Cum 



M 



M 

Z v> (>)Z v * r * (' " r ' )Z v >/ (' - r 2 ) 

_7=1 *=1 



AY 

s 

/=1 



Z V - r m(^^ 3 ) 



(114) 



Then, by the statistical independence of the source signals of assumption Al and 
repeatedly applying Cumulant Property 3, equation (114) reduces to 

C: (ti,r a> ii) =£c M /n|v,r.(/)v;r;(/-r l )v;r > (/-r 2 )v g .r,(/-r 3 ) J. (115) 

Using Cumulant Property 1 , the complex weights may then be pulled out in 
front of the cumulant operator in equation (1 15) to give 

< (T„r a ,r,) =EvXv;.v,C«/»[o(/)r;(/-r l )o(f-r J )r i (/-r 3 )]. (116) 

Arc y=l 

However, v v* = a 2 ri and equation (116) reduces 



- — i 

Cf (r„r 2 ,r 3 ) = £a^Cum[r, Wo (t-^Vj (t-t^t-z,)]. (117) 

Jrc 7=1 

From (1 17) it can be seen that spatial fourth-order cumulant matrix 3 in 
general can not be factored into Hermitian form, as was the case for spatial fourth- 
order cumulant matrix 1 and the spatial covariance matrix. However, if the 
"modified" steering vector elements are again defined as 

(118) 



^7 a rj V rj 



it can be factored in to bilinear form. 



cf(T 1 ,r 2 ,r 3 ) = V*Cj(r ls r 2 ,r 3 )V' 



(119) 



where the element in the r' h row and c' h column of the NxM "modified" mixing 
matrix V is 



[V] rc =v„. (1.2) 



Expanding equation (119), it is found that spatial fourth-order cumulant matrix 3 can 
be written as a sum of the outer products of the conjugate of the "modified" steering 
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V 

vector, v . , and steering vector scaled by the individual source signal's fourth-order 
cumulant. 

< fc^2.ti) = £c;(r |f r 2f T,)v>; (120) 

As before, it is yet to be proven if the mixing matrix V having full column 
5 rank is sufficient to guarantee that the "modified" mixing matrix V will have full 
column rank. However, it shall be assumed that the Hadamard product preserves the 
rank of the mixing matrix and therefore that the mixing matrix having full column 
rank is sufficient to guarantee that the "modified" mixing matrix has full column rank. 

If the "modified" mixing matrix has full column rank, by inspection of 
10 equation (120) it is clear that spatial fourth-order cumulant matrix 3 lies in the signal 
subspace spanned by the set of conjugated "modified" steering vectors. Again, like 
spatial fourth-order cumulant matrix 2, spatial fourth-order cumulant matrix 3 has no 
noise subspace. Note that in N. Yuen and B. Friedlander, "Asymptotic Performance 
Analysis of ESPRIT, Higher-order ESPRIT, and Virtual ESPRIT Algorithms," IEEE 
15 Trans. Signal Processing, Vol. 44, No. 10, Oct. 1996, pp. 2537-2550, as in H.H. 
Chiang and C.L. Nikias, "The ESPRIT Algorithm with Higher-Order Statistics," 
Proc. Workshop on Higher-Order Spectral Analysis, Vail, CO., Jun. 1989, pp. 163- 
168 and C.L. Nikias and A.P. Petropulu, Higher-Order Spectra Analysis: A Non- 
Linear Signal Processing Framework. (PTR Prentice-Hall, Upper Saddle River, NJ: 
20 1993), the elements/sensors are assumed to be omni-directional with unit gain so that 
aj=l. 

Spatial Fourth-Order Cumulant Matrix 3 Properties 

As with spatial fourth-order cumulant matrix 2, if the "modified" mixing matrix V 
has full column rank, spatial fourth-order cumulant matrix 3 will have many 
25 properties in common with spatial fourth-order cumulant matrix 1 . Properties 
associated with the development of a matrix-pencil and the associated separation 
technique with the assumption that the "modified" mixing matrix has full column rank 
are derived below. 

Property 1 : Spatial fourth-order cumulant matrix 3 is in general non- 
30 Hermitian. It will be Hermitian if and only if r, = r 3 = r and r 2 = 0 , i.e. C* (r, 0, r) 
and the sensors all have the identical gain for a given signal. 
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Property 2 : The trace of spatial fourth-order cumulant matrix 3 equals the sum 
of the signal fourth-order cumulants scaled by the sum of the sensor magnitudes to the 
fourth power. 

( n \ M N 

tr\C 4 t (r„r 2 ,r 3 ) -££aJcJ(r,,T 2 ,r s ) (121) 

V J Jm\ r =\ 

Property 3 : The column space of spatial fourth-order cumulant matrix 3, 



denoted as C 



(cf(r„r 2 ,r,)^ 



, is spanned by the set of conjugated "modified" 



steering vectors. 



sp 



'c[c*(T i9 T 29 T^ (122) 



Further, if the "modified" mixing matrix has full column rank, then the set of 
1 0 conjugated "modified" steering vectors are linearly independent and they form a basis 
for the column space of spatial fourth-order cumulant matrix 3. 

Property 4 : The rank of spatial fourth-order cumulant matrix 3 equals the rank 
of the mixing matrix, if V and V have full column rank. That is 

/>(c 4 ;(r 1 ,r 2 ,r 3 )j = >9(V) (123) 

15 if p{y) = p{v) = M > where /?( ) denotes rank. 

Property 5 : The "right" nullspace of spatial fourth-order cumulant matrix 3 
and the conjugate of the "left" nullspace of the mixing matrix are equal if the mixing 
matrix and "modified" mixing matrix have full column rank. 

N r [ci (r„r 2 ,r 3 )) = N (V) (124) 

20 The three definitions of the spatial fourth-order cumulant matrix presented all 

have the key attribute that they have no noise subspace. This feature allows us to 
avoid either having to use degrees of freedom in the array to estimate the noise 
subspace or having to make an assumption of temporal whiteness so the noise 
subspace will be gone at non-zero time lags. However, there are two main differences 

25 between definition 1 and definitions 2 and 3. 

First, definitions 2 and 3 have a computational advantage over definition 1 . 
This can be seen by comparing equation (76) with equations (93) and (111) where 
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definitions 2 and 3 require N 2 cumulants to be estimated while definition 1 requires 
Af 3 cumulants to be estimated. Second, while it remains to be rigorously proven that 
spatial fourth-order cumulant matrix 1 will have a rank equal to the number of signals 
if the mixing matrix has full column rank, it has not been proven that spatial fourth- 

5 order cumulant matrices 2 and 3 will have a rank equal to the number of signals if the 
mixing matrix has full column rank. This second difference arises from the fact that 
no proof that the Hadamard product preserves rank has been found. Therefore, it is 
assumed that for the special case of the "modified" mixing matrix it does and thus 
spatial fourth-order cumulant matrices 2 and 3 possess the derived set of properties 

10 needed to perform the blind source separation. However, if this assumption turns out 
not to be true, then unless the sensors in the array have identical manifolds, spatial 
fourth-order cumulant matrix 2 or 3 may not possess enough degrees of freedom to 
perform the separation. 

When the sensors in the array all have identical manifolds, the magnitude of 
1 5 the sensors' response, Iv^ = a i} , is a constant for each signal, that is 



= <Xij=<ij- (125) 



From equation (26) 

2>;=i ( 126 ) 

which, when all the sensors have identical manifolds, says that (125) can be 
20 substituted in to equation (126) to get 

2>j=M,j==l 027) 

and thus 

a,=j~. (.28) 
Further, if the manifolds are identical then 
25 0j= ^N (129) 
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for all j . Thus, for spatial fourth-order cumulant matrix 2 it was found that when the 
sensors have identical manifolds, substituting (129) in to equation (100) and factoring 
in to bilinear form leads to 

C:'(r l ,r 2 ,r } ) = i-VC:(r„r 2 ,r 3 )V H 

j (130) 

5 and therefore spatial fourth-order cumulant matrix 2 and spatial fourth-order cumulant 
matrix 1 are equivalent within a real scale factor. Following the identical path for 
spatial fourth-order cumulant matrix 3, one can find that 

Cf (r P r 2 ,r 3 ) = lv^(r 15 r 2 ,r 3 )V T 

" (131) 

= — cT(r lf * 2> r,) 

and therefore when the sensors have identical manifolds, spatial fourth-order 
10 cumulant matrix 3 is equivalent to the conjugate of spatial fourth-order cumulant 
matrix 1 within a real scale factor. Unfortunately, the assumption that all sensors 
have an identical spatial response is physically unrealizable and proves to be the 
undoing of the ESPRIT algorithm and its higher-order counterparts. 

Finally, in N. Yuen and B. Friedlander, "Asymptotic Performance Analysis of 
1 5 ESPRIT, Higher-order ESPRIT, and Virtual ESPRIT Algorithms," IEEE Trans. 

Signal Processing, Vol. 44, No. 10, Oct. 1996, pp. 2537-2550, it is claimed that the 
advantage of definition 3 over definition 2 is that when finite length data is used to 
estimate the spatial fourth-order cumulant matrix, spatial fourth-order cumulant 
matrix 3 will retain its Hermitian symmetry, subject to the conditions stated in 
20 property 1 . This property is not known to be of value to the matrix-pencil approach 
presented in the subsequent chapters and thus has not been evaluated for its validity. 

Provided below are a spatial fourth-order cumulant matrix pencil definition 
and associated properties. The indeterminacy of the standard eigen analysis of a 
spatial fourth-order cumulant matrix is described to motivate the use of a spatial 
25 fourth-order cumulant matrix-pencil by showing. The definition, properties, and 

spectral theory of matrix-pencils are presented including treatment of singular pencils 
and development of the novel concept of wide sense equivalence. The spatial fourth- 
order cumulant matrix-pencil is defined and its properties are derived for all three 
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spatial fourth-order cumulant matrix definitions. Finally, it is shown that the spectral 
analysis of the spatial fourth-order cumulant matrix-pencil provides a set of 
generalized eigenvectors that can be used as a basis for blindly separating the 
individual source signals. 

5 The formulation of a separation matrix to perform blind source includes 

finding a set of vectors that are each uniquely orthogonal to all steering vectors except 
one. The set of these vectors, possibly scaled by a normalization factor, form the 
columns of the separation matrix W that will diagonalize the mixing matrix V . The 
concept of blind source separation was previously described above, and the 

1 0 development of a technique to find a separation matrix using a spatial fourth-order 
cumulant matrix is provided below. 

A spatial fourth-order cumulant matrix signal subspace based spectral 
estimation technique is sought to perform the blind separation. In the context of a 
matrix subspace, spectral estimation implies eigen analysis and as such the terms 

1 5 spectral and eigen will be used interchangeably. In mathematics, eigenvalues are also 
often referred to as proper values. See, for example, P. R. Halmos, Finite- 
Dimensional Vector Spaces. (Springer- Verlag, New York, NY: 1987), which is 
hereby incorporated by reference in its entirety as if presented herein. Unfortunately, 
in general the standard spectral decomposition of the spatial fourth-order cumulant 

20 matrix will not directly provide a set of eigenvectors that will diagonalize the mixing 
matrix. A set of eigenvectors for the spatial fourth-order cumulant matrix does exist 
that will diagonalize the mixing matrix, however the entire signal subspace would 
have to be searched to find this unique set of eigenvectors. Thus, the standard eigen 
decomposition of the spatial fourth-order cumulant matrix possesses an indeterminacy 

25 that makes it undesirable. The indeterminacy of the standard eigen analysis is 
discussed in Section 5.2 in detail. 

The indeterminacy of the standard eigen decomposition can be overcome by 
the use of the generalized eigen analysis of a spatial fourth-order cumulant matrix- 
pencil. Therefore, a spatial fourth cumulant matrix-pencil of two spatial fourth-order 
30 cumulant matrices is defined by using two different sets of time lags, (0, 0,0) and 

( r P r 2> r 3)- 

Indeterminacy of the Standard Eigen Analysis : The formulation of a signal 
subspace based spectral estimation method to blindly separate the mixture of 
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statistically independent source signals is begun by examining the indeterminacy that 
exists in trying to use the standard eigenvectors of spatial fourth-order cumulant 
matrix 1 , C* (r, , r 2 , r 3 ) , to perform the blind source separation. In general, the results 

presented for spatial fourth-order cumulant matrix 1 are directly applicable to spatial 
fourth-order cumulant matrices 2 and 3 and thus will be presented for spatial fourth- 
order cumulant matrix 1 only. However, any differences or exceptions that may exist 
for spatial fourth-order cumulant matrices 2 and 3 will be appropriately noted. 

The standard eigenvalue problem for spatial fourth-order cumulant matrix 1 is 
defined as 

C 4 x (r p r 2 ,r 3 )e = >le. (132) 

A scalar X is said to be an eigenvalue of spatial fourth-order cumulant matrix 1 if the 
equality in equation (132) holds a non-null eigenvector e , an associated eigenvector. 
Rewriting equation (132), it is found that an eigenvector of spatial fourth-order 
cumulant matrix 1 lies in the "right" null space of a matrix-pencil. That is 

15 (Ct(h,T2>h)-M»)* = <> ( 133 > 

therefore 

eeN r (C 4 x (T x ,T 2 ,T 3 )-XI N ). (134) 

The matrix-pencil 

{C:(r 1 ,r 2 ,r 3 ),l}-c:(r 1 ,r 2 ,r 3 )->ll ;v (135) 

20 is non-singular with rank N even though C* (r p r 2 ,r 3 ) has rank M when X is not an 

eigenvalue. Thus, the eigenvalue X will reduce the rank of the matrix-pencil 
{C* (r„r 2 ,r 3 ),I w ) by a value rj &om , which is referred to as the "geometric 

multiplicity" of the eigenvalue X and is given by 

rj geom ^N-p(Cl{r XJ r 2y T 2 yXI N )^ (136) 

25 Since X is an eigenvalue if and only if the rank of the matrix-pencil is reduced, 

the eigenvalues can be found by searching for roots of the characteristic equation. 
That is, the values that cause the determinant of the matrix-pencil to be zero, 

det(c 4 x (r p r 2 ,r 3 )-^) = 0 (137) 

are eigenvalues. 
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The determinant in equation (137) is defined as the sum of products 
drtfC-fc.r,,*,)-*,^^ 

([C:(r„r 2 ,r 3 )]^-^A) 



(138) 



where the set {c } ,c 2 ,...,c N } is the /'* permutation of the first N positive integers, 

with the summation occurring over all L = N\ such permutations, and the scalar 
^represents the element in the r* row and c' h column of the identity matrix l N . 

The exponent in (138) is a scalar defined as a function of the set {c,,c 2 ,...,c w } by 



AM 



<l>,(c x ,C 2 ,...,C N ) = YjZn 



(139) 



where is the number of integers in the sequence c n+l ,...,c N that are less than c„ . 
Substituting equation (84) in to (138) 

( M 



det(C 4 x (r 1 ,r 2 ,r 3 )-^) = j;(-l) 



f M 



^7=1 

7=1 



(140) 



it becomes clear that each non-zero value of A that will cause the determinant in (140) 
to be equal to zero will be a linear combination of the individual source signal fourth- 
order cumulants. Therefore, since each eigenvalue is a linear combination of the 
source signals fourth-order cumulants it is reasonable to expect that the eigenvectors 
will be linear combinations of the associated steering vectors. This in fact can be 
shown to be true. Substituting equation (87) for the spatial fourth-order cumulant 
matrix, equation (132) becomes 

Z<(r„r 2 ,T 3 )v,v»e = Ae. (141) 



>1 



The vector inner product v**e will equal a scalar that is defined here as 



(142) 



then equation (141) becomes 



M £ 

e = 2-f c J( r i' r 2' r 3) v 7 



(143) 



53 



WO 03/090127 



PCTAJS03/11473 



Alternatively, it can be show that each eigenvalue is a linear combination of 

the source signal fourth-order cumulants by pre-multiplying (141) by the Hermitian 
transpose of the eigenvector and dividing by the inner product, e H e , since there is 

only one eigenvalue such that the equality in equation (141) for a particular 
5 eigenvector (See, for example, D. A. Harville, Matrix Algebra from a Statistician 's 
Perspective. (Springer-Verlag, New York, NY: 1999), which is hereby incorporated 
by reference in its entirety as if presented herein), each eigenvalue can therefore be 
expressed as 

v e H v v H e 

A-Zcjfe.^r,)— (144) 

7=1 e C 

10 which is clearly a linear combination of the individual source signals' fourth-order 
cumulants. 

If it is assumed that the mixing matrix V has full column rank, spatial fourth- 
order cumulant matrix 1 will have rank M by property 4. Therefore, there will be 
M non-zero eigenvalues whose sum will equal the trace of the spatial fourth-order 
1 5 cumulant matrix. That is, employing property 2 of spatial fourth-order cumulant 
matrix 1 

M M 

£^ =/r[c:(r I> T 2 ,r 3 )] = £c;(r 1 ,r 2 ,r 3 ). (145) 

Substituting (144) into equation (145) results in 

Z<(r l ,r 2 ,r 3 )S ! ^ i - = Z<(r„r 2 ,r 3 ). (146) 

20 Clearly, there exists an indeterminacy in the standard eigen analysis of spatial 

fourth-order cumulant matrix 1 . The same indeterminacy exists for definitions two 
and three and arises from the fact that in general, the identity matrix is not "similar" to 
spatial fourth-order cumulant matrix 1 . Therefore, a new matrix must be found, 
which for now can be referred to as the new B matrix, that is "similar" to the spatial 

25 fourth-order cumulant matrix to replace the identity matrix in the matrix-pencil and 
thus move to a generalized eigen analysis of spatial fourth-order cumulant matrix 1 
with respect to this new B matrix. Here "similar" means that the B matrix can be 
factored into a bilinear form, as was done for spatial fourth-order cumulant matrix 1, 
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« 



with the mixing matrix, and the modified mixing matrix for definitions 2 and 3, being 
two of the three factors and some diagonal matrix being the third. That is 

B = VDV H (147) 

where D is a diagonal matrix. 

5 Definition, Properties, and Spectral Theory of Matrix-Pencils 

A matrix-pencil is a special case of what is known in mathematics as a polynomial 
operator pencil. Polynomial operator pencils and their associated spectral problems 
arise naturally in many diverse areas such as differential equations, boundary value 
problems, control theory, harmonic systems analysis, wave propagation, elasticity 
10 theory, circuit simulation and modeling, and hydromechanics. See, for example, A. S. 
Markus, Introduction to the Spectral Theory of Polynomial Operator Pencils, 
Translation of Mathematical Monographs, Vol 71 . (American Mathematical Society, 
Providence, RI: 1988), which is hereby incorporated by reference in its entirety as if 
presented herein. In general, an n lh order polynomial operator pencil takes the form 

15 A(A) = A 0 +/IA, +••• + X* A„ (148) 

where X is a spectral parameter and the A/s are linear operators acting in a Hilbert 
space. A matrix-pencil, P (X) , is a first-order polynomial operator pencil that has the 

form 

V(X) = A-XB. (149) 

20 In general, matrix-pencils are classified as being either regular or singular. 

See, for example, A. S. Markus, Introduction to the Spectral Theory of Polynomial 
Operator Pencils, Translation of Mathematical Monographs, Vol. 71. (American 
Mathematical Society, Providence, RI: 1988), Z. Bai, J. Demmel, J. Dongarra, A. 
Ruhe, and H. van der Vorst, Templates for the Solution of Algebraic Eigenvalue 

25 Problems: A Practical Guide. (SIAM, Philadelphia, PA: 2000), K. Kanatani, 

Statistical Optimization for Geometric Computation: Theory and Practice. (Elsevier 
Science B.V., Amsterdam, The Netherlands: 1996), G. H. Golub and C. F. Van Loan, 
Matrix Computations. (The Johns Hopkins University Press, Baltimore, MD: 1996), 
F. R. Gantmacher, The Theory of Matrices, Volume I. (AMS Chelsea Publishing, 

30 Providence, RI, 1977), and F. R. Gantmacher, The Theory of Matrices, Volume II. 
(AMS Chelsea Publishing, Providence, RI, 1989), each of which is hereby 
incorporated by reference in its entirety as if presented herein. If the two matrices, 
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A, B , are square and the determinant of the matrix-pencil is not identically zero for all 
values of A , that is 

det (p(/l)) = det(A-AB)*0 (150) 

then the pencil is regular otherwise it is singular. Regular pencils have well defined 
5 eigenvalues which change continuously as functions of matrices A and B . Singular 
pencils, on the other hand, have eigenvalues that are discontinuous functions of A and 

B . Both types of pencils arise in practice and, are applicable to the BSS technique in 

accordance with the present invention. Note that the standard eigen problem is a 
regular matrix-pencil with B = I N . 

10 Properties and Spectral Theory of Regular Matrix-Pencils 

Regular pencils can be further classified as being Hermitian or non-Hermitian. Non- 
Hermitian matrix-pencils and their associated generalized non-Hermitian eigen 
problems arise when either A or B is non-Hermitian or B is not positive definite. 
Due to property 1 of the spatial fourth-order cumulant matrix, the spatial fourth-order 

15 cumulant matrix-pencil will in general be non-Hermitian. This will be shown to hold 
for all three definitions. Therefore, the focus will be on non-Hermitian regular pencils 
only and hence forth when the term regular pencil is used a non-Hermitian pencil is 
implied. Please see, for example, Z. Bai, J. Demmel, J. Dongarra, A. Ruhe, and H. 
van der Vorst, Templates for the Solution of Algebraic Eigenvalue Problems: A 

20 Practical Guide. (SIAM, Philadelphia, PA: 2000), K. Kanatani, Statistical 

Optimization for Geometric Computation: Theory and Practice. (Elsevier Science 
B.V., Amsterdam, The Netherlands: 1996), G. H. Golub, and C. F. Van Loan, Matrix 
Computations. (The Johns Hopkins University Press, Baltimore, MD: 1996), and F. 
R. Gantmacher, The Theory of Matrices, Volume I. (AMS Chelsea Publishing, 

25 Providence, RI, 1977), each of which is hereby incorporated by reference in its 
entirety as if presented herein, for discussions on Hermitian pencils. 

The characteristic polynomial of a regular N by N matrix-pencil, 

p (X) = det (P ( X)) = det ( A - /IB) (151) 

is by definition not identically zero for all values of X . The degree of p(X) is at 

30 most N . This means that there are N eigenvalues, which may be finite or infinite, 
with the roots of p (X) = 0 being the finite eigenvalues of the matrix-pencil. The set 

of eigenvalues of the matrix-pencil, P(^), are more commonly known as the 
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"generalized" eigenvalues of the matrix A with respect to the matrix B and are 
defined by 

A(A,B) = {*eD :det(A-zB) = 0}. (152) 

The eigenvalues of a regular pencil are continuous functions of A and B thus 

5 small changes in A and B cause small changes in the eigenvalues. If the degree of 

the characteristic polynomial is less than N , then the pencil is said to have 

N-M infinite eigenvalues, where M is the degree of the characteristic polynomial 
p(A) . The set of all eigenvalues, A(A,B), of a matrix-pencil is referred to as its 

spectra. See for example , and F. R. Gantmacher, The Theory of Matrices, Volume II. 

10 (AMS Chelsea Publishing, Providence, RI, 1989), which is hereby incorporated by 
reference in its entirety as if presented herein, and A. S. Markus, Introduction to the 
Spectral Theory of Polynomial Operator Pencils, Translation of Mathematical 
Monographs, Vol 71. (American Mathematical Society, Providence, RI: 1988), Z. 
Bai, J. Demmel, J. Dongarra, A. Ruhe, and H. van der Vorst, Templates for the 

1 5 Solution of Algebraic Eigenvalue Problems: A Practical Guide. (SIAM, Philadelphia, 

PA: 2000), K. Kanatani, Statistical Optimization for Geometric Computation: Theory 

and Practice. (Elsevier Science B.V., Amsterdam, The Netherlands: 1996), G. H. 

Golub and C. F. Van Loan, Matrix Computations. (The Johns Hopkins University 

Press, Baltimore, MD: 1996), and F. R. Gantmacher, The Theory of Matrices, Volume 

20 /. (AMS Chelsea Publishing, Providence, RI, 1977)[45-50], Note that as with the 
standard eigenvalues, an eigenvalue reduces the rank of the pencil by a value rf eom , 

which is referred to as the "geometric multiplicity" of the eigenvalue X . 

For each finite eigenvalue, any non-zero vector lying in the right null space of 
the matrix pencil evaluated at the eigenvalue is defined as a "right" eigenvector for 
25 that eigenvalue. 

eeN ,(A-/IB) (153) 
That is, for X e A( A,B) , any vector e that satisfies 

(A-AB)e = 0, e*0 (154) 

is an eigenvector corresponding to that eigenvalue. As with the eigenvalues of a 
30 matrix-pencil, the eigenvectors are often referred to as "generalized" eigenvectors. 
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For an infinite eigenvalue, any non-zero vector lying in the right null space of the 
matrix B is an eigenvector. That is any non-zero vector that satisfies 

Be = 0 (155) 

corresponds to an eigenvalue X = oo . An N by N regular matrix-pencil may not 

have N linearly independent eigenvectors. However, at least one independent 

eigenvector will exist for each distinct eigenvalue. As with standard eigenvectors, the 
set of generalized eigenvalues /l(A,B) is unique while the set of eigenvectors is not. 

Each regular matrix-pencil has two associated subspaces, denoted X and Y , 
that have the same dimension and satisfy 

Axer,Bx€7 VxeJf. (156) 

These subspaces are called right and left deflating subspaces, respectively. Further, 

span xeA , {Ax,Bx} = Y (157) 

and therefore 

AX +BX =Y . (158) 

1 5 Deflating subspaces are important in the development of techniques for solving 
regular generalized eigen problems (See, for example, P. Van Dooren, "Reducing 
Subspaces: Definitions, Properties, and Algorithms," Matrix Pencils, Proc. Pite 
Havsbad, Lecture Notes in Mathematics 973, Springer- Verlag, New York, NY, 1982, 
pp. 58-73, which is hereby incorporated by reference in its entirety as if presented 

20 herein)[5 1], such as the QZ algorithm, which is currently considered the most 
powerful technique for solving dense non-Hermitian regular generalized eigen 
problems. 

Finally, let X and Y be non-singular matrices where 

A = Y H AX B = Y H BX . (159) 

25 Then the matrix-pencil 

P(A) = A-AB (160) 

is "equivalent" to the matrix-pencil P(A) and X and Y are referred to as 
"equivalence transformations". The matrix-pencil P(A) has the same eigenvalues as 
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P(A) and its right eigenvectors, denoted e , are related to the right eigenvectors of 
the pencil P(/l) by the transformation 



e^X^e. (161) 
Properties and Spectral Theory of Singular Matrix-Pencils 

5 Matrix-pencils are singular if either they are rectangular or they are square and 

det (P(A)) = det(A-AB) = 0 VA (162) 

holds. For a square pencil to be singular, both matrices A and B must be singular and 
have a common nullspace. That is 

det(A) = det(B) = 0 

v ; w (163) 
N r (A)=:N r (B) 

10 are necessary and sufficient conditions for the pencil P(>l) to be singular. Both cases 

arise in practice and are significantly more difficult than regular pencils to handle. 
Since by definition the spatial fourth-order cumulant matrix is square, only square 
singular pencils are considered. 

The characteristic polynomial of a singular NbyN matrix-pencil is equal to 
15 zero for all functions of A . Therefore, the eigenvalues of a singular pencil are 

discontinuous functions of the matrices A and B and care must be taken in defining 
the eigenvalues of singular pencils. Obviously, the eigenvalues can no longer be 
found by finding the roots of the characteristic polynomial. Instead, the concept of 
reducing subspaces is used to define the eigenvalues and eigenvectors of a singular 
20 pencil. 

A pair of right and left reducing subspaces of the pencil P(>l) > denoted X 
and Y respectively, satisfy 

Ax€r,Bxer VxeJf (164) 

and 

25 span xe ^ {Ax,Bx} = Y . (165) 

where the dimension of the right reducing subspace, X , is greater than the dimension 
of the left reducing subspace Y by the dimension of the right null space of the pencil 
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over the field of all rational functions of X . That is 

dim (X ) = dim ( Y ) + dim (N r ( A - XB)) . (1 66) 

The reducing subspaces play a role similar to that of deflating subspaces for regular 
pencils. 

5 The rank of a matrix-pencil is in general a constant M for most values of X . 

p(p(X)) = M (167) 

However, for a certain set of values the rank of the pencil is "reduced", thus 
motivating the concept of a reducing subspace. The set of values, denoted /l(A,B) , 

that contain the values of X that cause the rank of the singular pencil to "reduce" are 
10 the eigenvalues or spectra of the singular pencil. 

l(A,B) = {z.e C:p(A-zB)<M} (168) 

The amount the rank of the pencil is reduced is the geometric multiplicity, rj^ om , of 
the particular eigenvalue. 

rj*™ =M-p(j>{X% Ae A(A,B). (169) 

1 5 Note that the eigenvalues of a singular matrix-pencil can be finite, infinite, or 
indeterminate. 

For each finite eigenvalue, any non-zero vector lying in the right null space of 
the matrix pencil evaluated at the eigenvalue is defined as a "right" eigenvector 
corresponding to that eigenvalue. 

20 eeiV r (A-/lB) (170) 

That is, for X e /l(A,B) , any vector e that satisfies 

(A-AB)e = 0, e*0 (171) 

is an eigenvector corresponding to that eigenvalue. For an indeterminate eigenvalue, 
any non-zero vector lying in the right null space of the matrix B , and thus in the right 
25 null space of matrix A , is an eigenvector corresponding to an indeterminate 
eigenvalue. Put in other words, any non-zero vector lying in the right nullspace 
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of A (or B ) is an eigenvector corresponding to an indefinite eigenvalue. Rewriting 
equation (171) and solving for X results in 

e H Ae 

*-7E- (,72 > 

Clearly, if e lies in the common right null space of A and B , X = 0/0 and thus the 
5 eigenvalue is indeterminate. 

As with regular pencils, the concept of "strict" equivalence is now defined. 
Let X and Y be non-singular matrices that do not depend on X where 

A = Y H AX B = Y H BX . (173) 

Then the matrix-pencil 

10 P(A) = A-^B (174) 

is "strictly equivalent" to the matrix-pencil P(/l) and X and Y are referred to as " 
strict equivalence transformations". The matrix-pencil P(>l) has the same 
eigenvalues as P(^). Also, the right and left reducing of the matrix-pencil P(/l), 
X and Y , are related to the right and left reducing subspaces of the pencil ¥(X) by 

15 X = X ] X y=Y H F. (175) 

The consequence of the effect of the equivalence transformation on the reducing 
subspaces is that the right eigenvectors of P (X) , denoted e , are related to the right 

eigenvectors of the pencil P(/l) by the transformation 

e = X"V (176) 

20 Wide Sense Equivalence of Matrix-Pencils 

The term "strict" equivalence was emphasized previously herein to distinguish it from 
what is defined here as "wide sense" equivalence. Given the M by N full row rank 
matrices X and Y that are independent of X , where 

A s Y H AX B e Y H BX (177) 

25 then the Af by N singular pencil P(>l)is said to be wide sense equivalent to the M 
by M non-singular pencil P {X) where 

P(A) = A-AB. (178) 
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Note that having rectangular X or Y is a sufficient condition to guarantee that P (/l) 
is singular. 

It will now be determined if the wide sense equivalence transformations 
preserve the spectra of P(A) and if the eigenvectors of P (X) are related to the 

5 eigenvectors of the non-singular pencil P(/l) by some equivalence transformation. 
To show this is true, let iei(A,B), i(A,B)c>l(A,B) and e beafiniteor 

infinite generalized eigenvalue and associated eigenvector of the singular pencil 

P (A) . Since the M by N matrix X has full row rank it has a "right" inverse. That 

is a N by M matrix X" ,r exists such that 
10 XX" lr = l M . (179) 

Also, since the matrix Y has full row rank, Y H has full column rank and therefore 
Y H has a "left" inverse. That is a M by N matrix (Y H )~ ! ' exists such that 

(Y h )~ 17 Y" =1^. (180) 

Clearly 

15 (Y H )""p(A)X- ,r =(Y ,I )""Y H (A-AB)XX- ,r =P(A). (181) 

The generalized eigenvalue problem is formulated as 

Ae = iBe. (182) 

Define the TV" by 1 vector y as 



25 



-[err 



(183) 



20 where e is an eigenvector of the matrix-pencil P(/l). The products y H Ae and 
y H Be are scalars, therefore 

i=4#- ( 184 > 

y Be 

Similarly, an eigenvalue of the non-singular pencil P(/l) is 



e H Ae 

^ = ^1^- (185) 
e"Be 
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Substituting (177) in to equation (184) results in 

. y H (Y H AX)e 

A = L Z1- » ( 186 ) 

y |Y BXje 

Then, using equation (1 83), (1 86) becomes 

. e H (Y H V ,/ (Y H AX)e e H AXe 
e H (Y H )" W (Y"BX)e 

5 Clearly every finite or infinite A e A ( A, B j is an eigenvalue of the pencil P (A) with a 
corresponding eigenvector 

e = X" lr e. (188) 

It can therefore be concluded that the set of finite and infinite eigenvalues of 
P ( A) , A ( A, B j is equal to the set of eigenvalues of the non-singular pencil. That is 

10 a(A,b) = A(A,B) (189) 

and that the eigenvector e of P (A) is related to the eigenvector of P(A) by the 
equivalence transformation 

e = X~ l 'e. (190) 

The Spatial Fourth-Order Cumulant Matrix-Pencil: Definition and Properties 

1 5 The spatial fourth-order cumulant matrix-pencil is defined on a pair on a pair of two 
spatial fourth-order cumulant matrices at delay lags (0,0,0) and (r,,r 2 ,r 3 ) , as 

P x (A,f)-C 4 x (0,0 5 0)-AC 4 x (r 1 ,r 2 ,r 3 ) (191) 

where the set of non-zero delay lags is denoted in vector form as f = [r, , r 2 , r 3 ] . As 

with the spatial fourth-order cumulant matrix, there are in three definitions for the 
20 spatial fourth-order cumulant matrix-pencil; each corresponding to the definition used 
for the pair spatial fourth-order cumulant matrices. Spatial fourth-order cumulant 
matrix-pencil 1 uses a pair of spatial fourth-order cumulant matrices 1 and is given in 
equation (191) above. Spatial fourth-order cumulant matrix-pencil 2 is defined on a 
pair of spatial fourth-order cumulant matrices using definition 2 as 

25 P; (A,f)s Cf (0, 0, 0) - AC*' (r, , r 2 , r 3 ) . (1 92) 
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Finally, spatial fourth-order cumulant matrix-pencil 3 is defined on a pair of spatial 
fourth-order cumulant matrices using definition 3 as 

P; (/L,f) = Cf (0, 0, 0) - AC*' (r x , r 2 , h ) . (193) 

Since it was shown in Chapter 4 that all three definitions have similar matrix 
5 properties if the Hadamard product preserves rank, the properties for spatial fourth- 
order cumulant matrix-pencil 1 are derived and any differences for spatial fourth- 
order cumulant matrix-pencils 2 and 3 noted. 

Spatial Fourth-Order Cumulant Matrix-Pencil Property 1 

Spatial fourth-order cumulant matrix-pencil 1 can be factored into Hermitian form as 
10 P x (A,t) = VP p (A,t)V h (194) 

where V is the mixing matrix and P r (A,f ) is an M by M matrix-pencil on a pair of 
diagonal signal cumulant matrices. 

P r (^3)^C:(0,0,0)-^(r 1 ,r 2 ,r 3 ). (195) 

Spatial fourth-order cumulant matrix-pencils 2 and 3 can be factored in to the 
1 5 bilinear form 

P;(A,f) = VP r (A,x)V H (196) 

and 

P;(A,T) = V-P r (A,x)V T (197) 

respectively. 

20 Spatial Fourth-Order Cumulant Matrix-Pencil Property 2 

The rank of spatial fourth-order cumulant matrix-pencil 1 equals the number of 
signals, M , for "most" values of X if V has full column rank. That is, for 
Ag^(C^(0,0 3 0),C^(r 1 ,r 2 ,r 3 )) and p(\) = M, 

p(P x (A,f)) = M. (198) 

25 Spatial fourth-order cumulant matrix-pencils 2 and 3 possess the same property if the 
Hadamard product (196) preserves rank. 
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Spatial Fourth-Order Cumulant Matrix-Pencil Property 3 

The spatial fourth-order cumulant matrix-pencil is in general non-Hermitian. It is a 
regular pencil if M = N and V has full column rank, otherwise, when M <N or if 
V does not have full column rank it is a singular pencil. Spatial fourth-order 
5 cumulant matrix-pencils 2 and 3 additionally require that the Hadamard product given 
in (4.65) preserve the rank of V for the pencil to be regular. 

Spatial Fourth-Order Cumulant Matrix-Pencil Property 4 

The spatial fourth-order cumulant matrix-pencil, P x (A,t)', is strictly equivalent to the 
regular pencil P r (>1,t) if P x (A,x) is a regular pencil. Otherwise, the spatial fourth- 
10 order cumulant matrix-pencil, P x (/I, t) , is wide sense equivalent to the regular pencil 
P r (>*,t) if the mixing matrix has full column rank. Spatial fourth-order cumulant 
matrix-pencils Pj[(A,t) and P*(A,f) additionally require that the Hadamard product 

given in equation (4.65) preserve rank. 

Spectral Analysis of the Spatial Fourth-Order Cumulant Matrix-Pencil 

1 5 The spectral theory of the spatial fourth-order cumulant matrix-pencil is explored in 
two manners. First exploiting the equivalence, strict or wide sense, it will be show 
that the finite spectra of the pencil P x (^,t) has a one-to-one mapping to set of signal 

fourth-order cumulants and thus each generalized eigenvalue can be associated with a 
source and its associated eigenvector with the signal's steering vector. The same 

20 relation is then shown by using the outer product expansion of the eigenvalue problem 
and the linear independence of the steering vectors. In both cases V is assumed to 
have full column rank and it is the finite generalized eigenvalues and their associated 
eigenvectors that are the focus since it is only the signal subspace that is of interest. 
As before, spectral theory for spatial fourth-order cumulant matrix-pencil 1 will be 

25 presented and any differences that may arise for spatial fourth-order cumulant matrix- 
pencils 2 and 3 will be noted. 

From property 4 of the spatial fourth-order cumulant matrix-pencil, P x (A, x) 
is equivalent, in a strict or wide sense, to the M by M regular pencil P r (A,t) . By 

the definition of equivalence, the set of finite and infinite spectra eigenvalues of 

30 
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15 



P x (A,f ) , i(C 4 (0,0,0),C 4 (*-„r 2 ,r 3 )) , is equal to the spectra of the pencil P r (A,t), 
A(c;(0,0,0),Cj(r„T 2 ,r 3 )). 

i(C 4 I (0 > 0,0),C 4 I (r 1 ,r 2 ,r 3 )) = A(C 4 r (0,0,0),C^(T 1 ,r 2 ,r3)) (199) 

m 

Clearly if P x (A,f ) is regular then 

i(C 4 I (0,0 5 0),C 4 x (r 1 ,r 2 ,r 3 )) = A(C 4 x (0 > 0,0),C 4 B (r l ,r 2 ,r,)). (200) 

Since the pencil P r is regular, its spectra can be determined by finding 
the roots of its determinate set equal to zero. Since P r (A,t) is diagonal by definition, 

the determinate will be the product of its diagonal components, i.e. 

det(P r (A,T)) = n«(0 > 0,0)-Ac 4 (T l5 r 2 ,r 3 )) 



(201) 



10 By inspection of (201 ), the spectra of P r (X, r) is the set 



A(C 4 (0,0,0),C 4 (r„T 2 ,r 3 ))=: 



c 4 (0,0,0) , 1 
z e €:: z = W e {l,. ..,M}\. (202) 



By Assumption A2 the signal cumulants are strictly non-zero and as such the spectra 

will contain M finite, non-zero eigenvalues, counting multiplicities, each 

corresponding to the ratio of a particular signal's cumulants at a zero and non-zero set 

of lags. Since the eigenvalues maybe complex, there is no fixed way to order them. 
For convenience, however, they will be ordered by their association with the j th 

signal as X } , i.e. 



(0,0,0) 

C rj ( T l> r 2> r 3) 



A; = 



(203) 



Of the M eigenvalues, there may only be K distinct values they assume, denoted as 
20 fi k , with each value having a multiplicity of rj k . Therefore, the eigenvalues can be 

grouped into K sets, denoted X k (C* (0,0,0),C 4 (r,,r 2 ,r 3 )), that contain the 

eigenvalues equal to ju k . 



X k (C 4 (0,0,0),C 4 (t„t 2 ,t 3 )) A (V C: ^ = A| 



(204) 
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Note that for a diagonal regular pencil, rj* om = tjf* = % . Clearly, 

A(C: (0,0,0),C: {t {M )) = A, (c: (0, 0, 0),Cj (r 1 ,r 2 ,r 3 ))u. .. 

u^(CJ(0,0,0),C:(r 1 ,r 2 ,r 3 )) 

For each /i A there are % signals with the identical ratio of their cumulant at zero lags 

to their cumulant at non-zero lags. This ratio is the inverse of what will be defined in 
5 the subsequent chapter as the normalized fourth-order auto-cumulant. 

A "right" eigenvector corresponding to a particular eigenvalue for a regular 
pencil will be a non-zero vector that lies in the right null space of the pencil evaluated 
at the eigenvalue. 

P r (A = >l y ,f)e 7 =0 (206) 

10 The eigenvector e y is a M by 1 vector with M-rj k zeros in the positions with 

indices corresponding to the columns of the diagonal pencil containing non-zero 
elements when evaluated at X = X } ,= fi k , a distinct eigenvalue. This is now illustrated 

with an example. 

For example, the M by M diagonal pencil, M > 3 , has the form 

c*(0,0,0)-/lc, 4 (r p r 2 ,r J ) 0 - - 0 



15 



0 



! C ;(0,0,0)-Ac;(r„r 2 ,r 3 ) i 

i \ 0 

0 - - 0 <(0,0 } 0)-Ac; ( (r„r 2 ,r 3 ) 



If 



< (0,0,0) _ < (0,0,0) 



Cj(*i.*2,*3) C r 3 ( T 1> T 2' T 3) 



then 



< (0,0,0) < (0,0,0) 

3) \- -77 r - "T7 \ ~ ^3 - Ml 
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and the distinct eigenvalue /i, has a multiplicity of two. The eigenvectors 
corresponding to the eigenvalues A, and Aj have the form 

e,=[e u 0 e 3J 0 o] T ;€{l,3} 

where e iy and e y are arbitrary scalars, both not equal to zero. Clearly the 

5 eigenvectors are not unique, however the locations of the non-zero elements are and 
the eigenvectors are uniquely associated with the signals. 

A set of M eigenvectors is required. However, when there are only K distinct 
eigenvalues, care must be taken in selecting the eigenvectors corresponding to a 
repeated proper value so that subsequent separation stages may be used to separate the 
10 signals associated with a repeated eigenvalue. The constraints for finding 

eigenvectors corresponding to repeated eigenvalues will be discussed in the following 
chapter. 

Consider now the spatial fourth-order cumulant matrix-pencil, P x (A,t) . If 

the mixing matrix V is assumed to have full column rank, then by spatial fourth- 
15 order cumulant matrix-pencil property 4, P x (A,t) is equivalent, either strictly or in a 

wide sense, to the pencil P r (A,t) . Thus, the pencil P x (A,t) has the same finite and 

infinite eigenvalues as P r (A,t) and the eigenvectors of P x (A,t) , designated e y , 

corresponding to the finite eigenvalues are related to the eigenvectors of P r (A,t) by 

e y =(v H )~\.. (207) 

20 Since the rows of V H are the conjugated signal steering vectors by definition, for 

V H (v H )" lr =I w (208) 

the columns of ( V H ) ~ 1r must be orthogonal to all steering vectors except the one in 

the row with the index corresponding to its own column index. Therefore, the 
columns of (V H ) are uniquely associated with a signal's steering vector. Thus it 

25 can be seen that 

ej V = { V H e,}" = {V H ( V" )"" e,}" = e, H . (209) 

Thus the eigenvectors, e y , can be used to blindly separate the signals. As 
mentioned earlier, the selection of the eigenvectors corresponding to a repeated 
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eigenvalue must be done with care. Clearly the output of e" V will be a linear 

combination of the rows of the mixing matrix and thus wilt result in a vector who is 

linear combination of the rows of the steering vectors of the signals associated with 
the eigenvector, that is the rows of V corresponding to the non-zero elements of e y . 

5 For spatial fourth-order cumulant matrix-pencils 2 and 3, if the Hadamard product 
preserves the rank of V , then a similar result is found by exploiting a corresponding 
equivalence transformation between the "left" eigenvectors of the pencils P r 

and P x (A, t) . Note that similar results will be achieved with the left eigenvectors, 
i.e. those 1 by N vectors, dj , that lie in the left null space of the pencil evaluated at 
10 the eigenvalue. 

The preceding spectral analysis of the spatial fourth-order cumulant matrix- 
pencil can be performed in an alternative way by expanding the pencil into a sum of 

outer products and exploiting the linear independence of the steering vectors. The 
generalized (right) eigenvalue problem for the matrix-pencil P x (^,t) is defined as 

15 C£ (0,0,0ft -ACj^r^,)^ (210) 

Substituting equation (87) for the spatial fourth-order cumulant matrix in (210) and 
rearranging terms results in 

£ v . k (0, 0, 0) - tx k c\ (r, , r 2 , r 3 )] v % = 0 . (211) 

7=1 

If the mixing matrix is assumed to have full column rank its columns are linearly 
20 independent and thus 

[e* (0,0,0)-^; (r p r 2 ,r 3 )]v;e 4 =0 (212) 
for all j for the equality in equation 21 1 to hold. This leads to the 

Mk =—, 7 ( 213 > 

C ry ( T l' r 2' r 3) 

for any j when v"e t *0. Since by property 5 V H and (r, , t 2 , r 3 ) have a 

25 common right null space, any eigenvector that lies in the right null space of V H has a 
corresponding indeterminate eigenvalue, since e k can only lie in the null space of V H 

if the pencil is singular. Therefore, as with the equivalence approach presented 
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previously, the eigenvalues and their associated eigenvectors are uniquely associated 
with the source signals with the ratio 

cj (0,0,0) 



*i(*i,r 2 ,r 3 ) 



(214) 



acting as a discriminant. Similar results hold for spatial fourth-order cumulant 
5 matrix-pencils two and three if the modified mixing matrix V has full column rank. 

A blind source separation technique in accordance with the present invention 
is described utilizing all three spatial fourth-order cumulant matrix definitions and the 
conditions of identifiability are discussed. A normalization method is also developed 
to allow the separation algorithm to maximize the separation power efficiency (SPE). 
1 0 The concept of a normalized fourth-order auto-cumulant is presented and a method of 
selecting eigenvectors for repeated eigenvalues to facilitate cascaded processing is 
developed. 

Figure 5 is a functional block diagram of a blind source technique and 
processor in accordance with an embodiment of the present invention. Blindly 
1 5 separating M statistically independent source signals under the narrowband 

assumption requires finding a N by M separation matrix W that will diagonalize 
the mixing matrix V . That is, from equation (39), a separation matrix is sought such 
that 



W H V = 



p x 0 0 

0 p 2 \ : 

1 \ \ 0 
.0 - 0 p M \ 



(215) 



20 Calculating this separation matrix is one function of the blind source separation 
technique. 

It was previously explained that the generalized eigenvectors of the spatial 
fourth-order cumulant matrix-pencil would separate signals based on the ratio of their 
cumulants at two different sets of time lags. These results are exploited in the 
25 formalization of a technique to calculate a separation. This technique will 

theoretically minimize the residual interference-to-signal ratio (ISR) and, with the 
proper normalization, to maximize the separation power efficiency (SPE). 
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As also previously described, a factor in finding the eigenvectors of the spatial 
fourth-order cumulant matrix-pencil is what was referred to as the normalized fourth- 
order auto-cumulant. This arises from the finite eigenvalues of the spatial fourth- 
order cumulant matrix-pencil equaling the inverses of the individual source signals 

5 normalized fourth-order auto-cumulants; the aforementioned ratio of cumulants at two 
different sets of time lags with the set of time lags in the numerator cumulant equal to 
zero. Specifically, since it was shown that if a signal has a unique normalized fourth- 
order auto-cumulant, the associated eigenvector is orthogonal to all of the other 
signals' steering vectors, the normalized fourth-order auto-cumulant can be thought of 

10 as acting as a signal discriminator function. In this chapter the normalized fourth- 
order auto-cumulant will be defined and some comments on its use as a discriminator 
function will be made. 

When multiple signals have the same normalized fourth-order auto-cumulant 
at the set of delay lags used, repeated eigenvalues occur. In order to facilitate 

15 repeating the separation technique at a new set of time lags where ideally the signals 
will no longer have equal normalized fourth-order auto-cumulants, care must be taken 
in selecting the eigenvectors of a repeated eigenvalue in order to guarantee the 
resultant steering vectors remain linearly independent. Criteria are presented below 
for selecting eigenvectors associated with a repeated eigenvalue. It will be shown that 

20 the resultant set of new steering vectors remains linearly independent. 

One measure of performance used in assessing the blind source separation 
algorithm is separation power efficiency (SPE). To achieve an SPE of 1, the inner 
product of the separation vector and the associated steering vector must have a 
magnitude of 1 . To achieve this, the eigenvector, which forms the basis of the 
25 separation vector, must be normalized since, although it is co-linear with the steering 
vector, it does not necessarily have the correct magnitude to guarantee that an SPE of 
1 can be achieved. Therefore, a normalization algorithm is developed to guarantee 
that the SPE is maximized. Since there are three definitions for the spatial fourth- 
order cumulant matrix, different normalization techniques will be required for each. 

30 Conditions of identifiability are presented below that will allow an 

instantaneous linear mixture of signals to be separated. These include the linear 
independence of the steering vectors, the statistical independence and non-Gaussian 
nature of the source signals, and existence of a set of time lags where every signal has 
a distinct normalized fourth-order auto-cumulant, just to name a few. 
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Finally, the spatial fourth-order cumulant matrix-pencil based algorithm is 
presented step by step in a flow diagram. Each step is discussed and critical issues, if 
any, are presented. Difference that may exist in using the different spatial fourth- 
order cumulant matrix definitions in the algorithm will be highlighted. 



5 Normalized Fourth-Order Auto-Cumulant as a Signal Discriminator 

It was previously shown that the generalized eigenvalues of the spatial fourth- 
order cumulant matrix-pencil are 

c 4 (0,0,0) 

*,« 4 y v ;'euu.,M. (216) 

For separation to occur, a distinct eigenvalue, Xj , is required for each signal. 
10 Therefore, Aj acts as a "signal discriminant". To investigate this discriminant and its 
properties, the normalized fourth-order auto-cumulant of the f h signal is defined as 

S fa*..'.)- c ; (0>0-0) • < 217 > 

Clearly the generalized eigenvalue associated with the signal is the inverse of its 
normalized fourth-order auto-cumulant. 

15 ^j=znrr~ r ( 218 > 

By assumption Al, the signal r j (t) is a stationary, non-Gaussian random 

process. Further, by assumption A2, the signal is assumed to have a zero mean with 
power Pj , and a non-zero fourth-order moment. These assumptions ensure the 

signals fourth-order cumulant exists and is non-zero. A necessary extension to these 
20 assumptions is that the set of time lags, (r, , t 2 , r 3 ) , is chosen so that the fourth-order 

auto-cumulant also exists and is non-zero. Thus it may be assumed that the 
normalized fourth-order auto-cumulant exists and is finite, 

Since the signal is assumed to be a stationary random process, its moments are 
dependent only on the time differences or lags (r,,r 2 ,r 3 ) . Therefore, the normalized 

25 fourth-order cumulant is a three dimensional function. Thus, for separating signals 
there are three independent variables to work with to ensure the signals have unique 
normalized fourth-order auto-cumulants as compared to a second-order technique that 
has only one independent variable. This is yet another distinct advantage of the 
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fourth-order cumulant approach over the second-order spatial correlation based 
approach. 

The normalized fourth-order auto-cumulant will, in general, be complex 
valued. Although the signal's cumulant at the set of time lags (0,0,0) will have a 

5 real value, the cumulant at time lags (f X9 r 29 r 3 ) will be complex valued. Thus, the 

normalized fourth-order auto-cumulant will contain phase information that is a 
function of (r, , r 2 , r 3 ) . Source signal emitter attributes such as phase noise, carrier 

frequency, transmit filter response, amplifier response, transmitter clock jitter, 
propagation channel transfer function, etc., will contribute to the source signal's 
10 normalized fourth-order auto-cumulant. From the definition of received source 
signal, 

oW-V^/M (219) 

it is clear, employing cumulant property 1 , that the normalized fourth-order auto- 
cumulant is not a function of the signal's power. 



15 



_4, ,_<0W*3) 



(220) 



Ctim[r,(f)r / t (r-r 1 ) 0 (/-r 2 )r;(/-r 3 )] 

_ PjCum [mj(t) m) (t - r, ) m j {t - 1 2 ) m) (t - r 3 )] 
PjCum[mj (t)m) (t) mj (t)m) (/)] 

= Cum [ntj (t ) m) {t - r, ) m j (t - r 2 ) m] (t - t 3 )] 
Ctm[mj(t)mi(t)mj{t)m]{t)\ 

Thus the signals are discriminated by having an underlying waveform that has a 
unique normalized fourth-order auto-cumulant, not by having different signal powers. 

As stated above, the unit power modulated signal from the /* emitter is 
affected by transmitter unique properties. Since in practice two transmitters hardly 
20 ever produce identical signals, a source signal will most likely have a unique fourth- 
order auto-cumulant function and therefore have a unique normalized fourth-order 
auto-cumulant function. Therefore, it is expected that a set of time lags will exist 
where a group of signals will have unique normalized fourth-order auto-cumulants 
and can therefore be separated. 
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Selecting Eigenvectors for Repeated Eigenvalues 

When multiple signals' normalized fourth-order auto-cumulant functions 
happen to have the same value at the set of time lags chosen, the problem of repeated 
eigenvalues arises. In this situation care must be taken in selecting the associated 
5 eigenvectors to ensure that the set of separation vectors formed from these 

eigenvectors will transform the mixing matrix into a new reduced dimension mixing 
matrix with full column rank. This guarantees that the separation algorithm can be 
repeated to separate the signals associated with the repeated eigenvalue that are mixed 
by this new reduced dimension mixing matrix. 

1 0 Figure 6 depicts repeating the separation algorithm for a single repeated 

eigenvalue. In Figure 6, the M by 1 vector output of the first separation stage, y (t) , 

is partitioned into two vectors, one with the M-rj K elements of y (t ) that correspond 

to unique eigenvalues, denoted as y s (f ) , and the other, denoted (f ) , with the rj K 

elements of y (f) that correspond to the repeated eigenvalue, fi K , which has 

15 multiplicity rj K . As in the first separation stage, a new tj K by tj k separation is 

sought that will separate the signals associated with the repeated eigenvalue. The new 

separation matrix is found by repeating the spatial fourth-order matrix-pencil 

algorithm at a different set of time lags. However, to use the algorithm a second time, 

the conditions for identifiability discussed in a following section must hold. One of 

20 these is that the mixing matrix has full column rank and therefore the new reduced 

dimension, rj K by 7] K mixing matrix must have full column rank. Each repeated 

eigenvalue would have a similar new mixing matrix and would have the separation 
algorithm repeated. 

The requirements for selecting the eigenvectors associated with a repeated 
25 eigenvalue so as to guarantee the new mixing matrix, , will have full column rank 

are derived below. The set of integers that represent the signals associated with a 
particular eigenvalue is defined as 

g k ={je{\A:.,M}:lj =//,}• (221) 

Recalling that the eigenvalue Xj is equal to the inverse of the f h signals normalized 
30 auto-cumulant. 

c r (0,0,0) 

V-r — h (222) 
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Since there are M signals there will be M eigenvalues, including repetitions, of 

which only K are distinct. The distinct eigenvalues are denoted as 

/u k9 k € 1, 2, . . . , K . The scaled or normalized M associated N by 1 eigenvectors 

are the columns of the separation matrix W . 

5 "j=rA (223) 

From equation (215) it can be seen that if the separation matrix separates all the 
signals, the resultant matrix product W H V is diagonal. This occurs only if there are 
M distinct eigenvalues. As previously described, for a repeated eigenvalue, the 
eigenvectors separated the signals associated with a particular eigenvalue from those 
10 not associated with it; however, the resultant scalar was a linear combination of the 
associated signals. Therefore, the matrix product W H V will have rj k rows, indexed 

by the set g k , each with rj k non-zero elements in the columns indexed by the set g k 
for k = 1, 2, . . . , K . This holds for the M distinct eigenvalues case as well in which 
case K = M and rj k = 1 , V£ , k = j , and each g k has only one element in it, k = j, 

15 and thus W H V is diagonal. 

The new mixing matrix V M consists of the 7j k columns of W H V indexed by 
the set g k with the rows indexed by the integers not in the set g k removed, i.e. the 
rows with all zeros. Thus is a rj k by rj k matrix that linearly mixes the rj k signals 
associated with the repeated eigenvalue to form the tj k elements of the vector (t) . 
20 Since is to the second stage separation process as V is to the initial separation 

process, it must have the same properties as V , primarily it must have full column 
rank. 

P(V^) = % (224) 
The condition for selecting w 7 , je g k that ensures has full column rank must 

25 now be determined. 

Since w y is a scaled version of the eigenvector e y , a constraint on Wj is in 
fact a constraint on the selection of e, . It will now be shown that requiring the 
eigenvectors , for j e g k to be linearly independent is sufficient to guarantee that 
will have full column rank if V has foil column rank. 

30 
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10 



Separation Vector Formation: Normalizing the Eigenvectors 

Although the eigenvectors are orthogonal to all signals except those with 
normalized auto-cumulants equal to the associated eigenvalue, the inner product 



(225) 



does not guarantee that the maximum SPE of 1 will be achieved. Therefore, the 

separation vectors are formed by scaling each eigenvector with a real normalization 
factor Yj that guarantees a maximum SPE of 1 . 



w . = r,e, 



(226) 



For repeated eigenvalues, the normalization factor will have a different effect and 
there is no clear advantage to normalizing the eigenvectors associated with repeated 
eigenvalues in the first separation stage. Further study of the effect of repeating the 
separation algorithm on the achievable SPE is needed. 



From equation (39) it can be seen that the inner product 



H 



(227) 



15 results in a "loss" term p, that in general is complex. From equation (63) it can be 



seen that the SPE is 



£=w, H v,v>,. 



(228) 



Inserting equation (227) in to (228) results in 



Zj = PjPj = \pj 



(229) 



20 For the SPE to be 1 requires 



Pj 



= 1 



(230) 



and thus a normalization factor y } is needed such that 



H _ . 

W • v ■ 
J J 



* H 1 

Yj*j y j\-y) 



(231) 



and therefore 



25 



1 



(232) 
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Calculation of this scale factor for the eigenvectors associated with unique 
eigenvalues is dependent on the particular spatial fourth-order cumulant matrix 
definition used since the available variables are different. For spatial fourth-order 
cumulant matrix 1 a normalization factor that will guarantee the SPE will achieve a 

5 maximum of 1 can be found and will be derived in the following sub-section. 

However, it will be shown for spatial fourth-order cumulant matrices 2 and 3 that the 
existence of the modified mixing matrix V in their bilinear forms causes the set of 
equations available for solving for the normalization factor to be under specified and 
thus using definitions 2 and 3 cause a power loss preventing the separation from 

1 0 achieving an SPE of 1 . 

Normalization for Spatial Fourth-Order Cumulant Matrix 1 

The only data available to work with in normalizing the eigenvectors are the 
spatial fourth-order cumulant matrices, the eigenvalues, and associated eigenvectors. 
From equation (225), 

15 £j = *"vj ( 233 > 

and as previously described, it is known, assuming is associated with a distinct 
eigenvalue, that 

V H e,=[o ... 0 €j 0 ... 0] T (234) 

with the non-zero element in the /* position. Therefore, since spatial fourth-order 
20 cumulant matrix 1 can be factored into Hermitian form as 

C 4 X (0,0,0) = VC* (0,0,0) V H (235) 

where C* (0,0,0) is an M by M diagonal matrix, the product 

C 4 x (0,0,0)e 7 =VC;(0,0,0)V H e y 



= VC;(0,0,0){[0 ... 0 e) 0 ... 0] T J 

= vj[0 ... 0 e)c\ (0,0,0) 0 ... 0] T J 
= e)c\ (0,0,0) v,. 



(236) 
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15 



The Euclidian or L norm of (236) is then 



|c 4 x (0,0,0)e y | 2 = |^(0,0 ) 0)v > 

^ (0>°,0) ||v, 



(237) 



However, since 



-4 
=i 



(238) 



then 



||C 4 ,(0,0,0)e ; || 2 =|^||c;(0,0,0)|. 



(239) 



Further, pre-multiplying the product in equation (236) by the Hermitian transpose of 
the /* eigenvector results in the scalar 



e, H c;|(0,0,0)e, =e/ rj (0,0,0)e, H v, 

= s j£j c A rj (0,0,0) 

< (0,0,0). 

Taking the ratio of the absolute value of (240) to (239) produces the scalar 

2 



(240) 



e y H C^(0,0,0)e, 



K (o,o,o)e y .| 



;(o>o>o)i _. 



(o>o,o)| 



(241) 



Thus the unknown denominator in (232) has been solved for and therefore the 
normalization factor, when using spatial fourth-order cumulant matrix 1 , is 

J_ lcI(o,o,o)i ; ||, 

f j ~ I _ I - 



eJC*(0,0,0)e, 



Normalization for Spatial Fourth-Order Cumulant Matrices 2 and 3 



(242) 



Spatial fourth-order cumulant matrix 2 and 3 are not factorable into Hermitian 
form but instead factor into a bi-linear form as shown previously as 



C 4 ,'(r 1 ,r 2 ,r 3 ) = VC:(r„r 2 ,r 3 )V H 



(243) 
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and 



Cf (r„r 2 ,r 3 ) = V*C; (r„r 2 ,r 3 ) V T (244) 



respectively. From the results in chapter 5 it is clear that pencils formed using spatial 
fourth-order cumulant matrices 2 and 3 will have the same eigenvalues with the 
5 associated eigenvectors of spatial fourth-order cumulant matrix-pencil 3 equal to the 
conjugate of the eigenvectors of spatial fourth-order cumulant matrix-pencil 2. Thus, 
since the normalization factor is real, if it exists it would be the same for both 
definitions. 

Since spatial fourth-order cumulant matrices 2 and 3 are not factorable into 
10 Hermitian form, the modified mixing matrix has to be dealt with in attempting to 
estimate the normalization factor given in equation (232). Unfortunately, in general 

ejfv, * Sj . (245) 

Further, even for distinct eigenvalues, the eigenvector, specifically the right 
eigenvector , is in general no longer orthogonal to all modified steering vectors 

1 5 except Vj and it is not guaranteed that the modified steering vectors have a Euclidian 
norm of 1 , that is in general 

vjv, * 1 . (246) 

Thus the properties that spatial fourth-order cumulant matrix 1 possessed that allowed 
for the estimation of \€j are not shared by spatial fourth-order cumulant matrices 2 

20 and 3. 

Solving for given only spatial fourth-order cumulant matrix 2 or 3, the 

generalized eigenvalues of the associated pencil, and their associated left, d 7 , and 
right, e , eigenvectors is under specified. Since spatial fourth-order cumulant matrix 

1 can be factored into Hermitian form, the left and right eigenvectors of spatial fourth- 
25 order cumulant matrix-pencil 1 are related by a Hermitian transformation, that is if 

is a right eigenvector of spatial fourth-order cumulant matrix-pencil 1 then e* 1 is a left 

eigenvector. Thus the product e"C* (0,0,0)^ has only two unknowns, one 



being Ej . Similarly, ||c* (0,0,0)e y || 2 has the same two unknowns since 



= 1 and 
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therefore s can be solved for. However, for spatial fourth-order cumulant matrix 



2, and similarly for 3, 



C 4 x >,0,0)e y = e]c A rj (0A0)\j 

2 

; (<W) Nl 



(247) 



and 



d,C 4 ; (0,0, 0)e,. = eft (0,0,0)d,v, 

= s)8/ rj (0,0,0) 

where d, is the 1 by N left eigenvector associated with the /* signal and 

Sj = djVj . 



(248) 



(249) 



Thus there are two equations and four unknowns. Attempting to establish four 
equations by using the spatial fourth-order cumulant matrix 2 at delay lags (t 19 t 2 ,t 3 ) 

1 0 results in now having four equations and five unknowns. Thus, for spatial fourth- 
order cumulant matrix 2, and similarly for 3, solving for the normalization factor 
results in an under specified problem that cannot be solve. This is yet another 
advantage of definition 1 , since a normalization factor that guarantees a maximum 
SPE of 1 can be solved for. 

15 Conditions for Identifiabilitv 

Identifiability deals with the ability of the blind source separation algorithm to 
uniquely associate a separation vector with an arriving source signal and thus separate 
it out of the linear mixture by suppressing the other signals. In order for the proposed 
blind source separation algorithm to perform the separation certain conditions must be 

20 met. Some have already been presented as signal and noise assumptions and are 
restated here as conditions that are imposed for the spatial fourth-order cumulant 
matrix-pencil based blind source separation algorithm to achieve the separation of the 
source signals. The fewer conditions for identifiability required the more powerful 
the algorithm will be in the sense of handling a wider variety of source separation 

25 problems. Five conditions for identifiability, CI1 through CIS, are given below, 
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CI1 : The mixing matrix, V , has full column rank. This requires the number of 
sources be less than or equal to the number of sensors, i.e. M <N , and that the signal 
steering vectors be linearly independent. 

CI2: The normalized fourth-order auto-cumulant, c * (r, , r 2 , r 3 ) , is different for each 

5 signal. The algorithm may be repeated at different sets of time lags (t 19 t 29 t 3 ) in a 

second separation stage operating only on the signals in the first stage that had 
identical normalized fourth-order auto-cumulants. 

CD: The M source signals illuminating the array are statistically independent non- 
Gaussian stationary random processes and are stationary to order four over the spatial 
1 0 fourth-order cumulant matrix estimation period. 

CI4: The noise processes are stationary Gaussian random processes. They need not 
be spatially or temporally white. Stationarity is required only over the estimation 
period of the spatial fourth-order cumulant matrix. 

CIS: For spatial fourth-order cumulant matrix definitions 2 and 3, the Hadamard 
1 5 product 

V=V0V0V (250) 

preserves the rank of the mixing matrix V , i.e. the modified mixing matrix has full 
column rank. This condition is not required when using spatial fourth-order cumulant 
matrix definition 1 . 

20 Separation Matrix Formation Algorithm 

Figure 7 and Figure 8 are flow diagrams of a process for performing the blind 
source separation using the spatial fourth-order cumulant matrix-pencil in accordance 
with an embodiment of the present invention. The algorithm requires as inputs the set 
of time lags (t { ,t 29 t 3 ) where (r,,r 2 ,r 3 )* (0,0,0). Delay lag values, ti,t 2 ,t 3 , are 

25 provided at step 61 and senor data values, x(t), are provided at step 63. It is 

recommended, that either r, ^r 2 or r 3 * 0 in order to preserve phase information. 

This will reduce the occurrence of repeated eigenvalues and thus reduce the number 
of times the separation has to be repeated. 

At step 60, the estimation of the spatial fourth-order cumulant matrix at lags 
30 (0,0,0) and (r,,r 2 ,r 3 ) is performed matrix element by matrix element. Since a 

cumulant cannot be estimated directly from the data, all moments up to order four 
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must be estimated. The estimation can be done either in a real time manner as the 
array samples the propagating wavefields and generates the sensor output data x(/) 

or after the entire data set required to estimate the cumulants has been captured. 

After the spatial fourth-order cumulant matrices C* (0,0,0) and C*(r„r 2 ,r 3 ) 
5 have been estimated, the generalized eigen analysis of the matrix-pencil P % is 
performed at step 62 to determine its finite spectra, X[C\ (0,0,0),C* (T,,r 2 ,r 3 )) . At 

step 64, the number, K, of distinct finite eigenvalues and multiplicity of each distinct 

eigenvalue are determined. The spectra will have M finite values, counting 

multiplicities, each corresponding to the normalized fourth-order auto-cumulant of a 
1 0 signal. Of the M eigenvalues there will be K distinct eigenvalues, fi k , 

k e 1,2,. . . , K , each with a multiplicity rj k . For each distinct eigenvalue rj k linearly 

independent eigenvectors are calculated. An index, k, is set to zero at step 66. The 

index, k, is used to ensure that each distinct eigenvalue is addressed. At step 68, the 

distinct eigenvalue index, k, is compared to the number of distinct finite eigenvalues, 

1 5 K. Given that at least one distinct eigenvalue exists, K will not be equal to zero. 

Thus, on the first iteration, k will be less than K, and the process will proceed to step 

72, as indicated by the circled letter "A" in Figure 7 and Figure 8. At step 72, it is 

determined if the multiplicity, T| k , is greater than 1 . If the multiplicity, r| k , is not 

greater than 1, the process proceeds to step 74. At step 74, an eigenvector ej, is 

20 calculated for the eigenvalue for the k* distinct eigenvalue (Xj = jik). For each 

JL = fx k with a multiplicity of one, a normalization factor y s is calculated at step 76. 

The separation vector is formed as w, = yjk j at step 78. The separation vector, Wj, is 

utilized (appended) to form the separation Matrix, W, at step 80, wherein the 
separation vectors are (by definition) the columns of the separation matrix, W. After 

25 the separation vector, Wj, is appended to the separation matrix, W, the index, k, is 
incremented at step 82. The process then proceeds to step 68, as indicated by the 
circled letter "B" in Figure 7 and Figure 8. At step 68, k is compared with K. If k is 
greater than K, then the separation matrix W is provided and available for subsequent 
processing at step 70. If k is not greater than K (step 68) then the process proceeds to 

30 step 72, as indicated by the circled letter "A" in Figure 7 and Figure 8. At step 72, it 

is determined if the multiplicity, rjk, is greater than 1 . If the multiplicity, is greater 

than 1 , the process proceeds to step 84. At step 84, linearly independent 

eigenvectors eigenvector ej, are calculated for the distinct eigenvalues (Xj = |ik). For 
each repeated eigenvalue the rj k separation vectors are set equal to its associated 
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eigenvectors as w y = e 7 , at step 86. The separation matrix, W, is formed at step 80 

by appending the separation vectors Wj. The index, k, is incremented again at step 
82, and the process is repeated until all distinct eigenvalues are addressed (k is greater 
than K at step 68). At step 68, if k is greater than K, the separation matrix, W, is 
5 provided and available for subsequent processing at step 70. At step 71, the 

separation matrix, W, is multiplied by the input signal x(t) to perform the separation. 
More specifically, a matrix representation of the input signal x(t) is multiplied by the 
Hermitian transpose of the separation matrix, W H , in accordance with the following 
equation. 

10 y(/) = W H x(r) (251) 

A BSS technique as described herein may be embodied in the form of 
computer-implemented processes and system for practicing those processes. A BSS 
technique as described herein may also be embodied in the form of computer program 
code embodied in tangible media, such as floppy diskettes, read only memories 

1 5 (ROMs), CD-ROMs, hard drives, high density disk, or any other computer-readable 
storage medium, wherein, when the computer program code is loaded into and 
executed by a computer, the computer becomes a system for practicing the invention. 
The BSS technique as described herein may also be embodied in the form of 
computer program code, for example, whether stored in a storage medium, loaded 

20 into and/or executed by a computer, or transmitted over some transmission medium, 
such as over the electrical wiring or cabling, through fiber optics, or via 
electromagnetic radiation, wherein, when the computer program code is loaded into 
and executed by a computer, the computer becomes a system for practicing the 
invention. When implemented on a general-purpose processor, the computer program 

25 code segments configure the processor to create specific logic circuits. 

A blind source separation (BSS) technique in accordance with the present 
invention provides a robust higher-order cumulant based principle component blind 
source separation technique that performs well at low signal-to-noise ratios with the 
potential of performing well in the presence of temporally and spatially correlated 
30 noise. Furthermore, a new definition of a spatial fourth-order cumulant matrix suited 
to blind source separation with non-equal gain, directional sensors is provided, the 
definition of a spatial fourth-order cumulant matrix-pencil using temporal information 
is provided, the concept of the separation power efficiency as a measure of the 
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algorithm's performance is provided, and the concept of wide sense equivalence 
between matrix-pencils are also provided. 

Applications of the BSS technique in accordance with the present invention 
include spectral monitoring, be it for signal intelligence or other applications such as 

5 radio astronomy, where Gaussian random noise processes dominate the received 
signals out of an array. This fourth-order array signal processing BSS technique in 
accordance with the present invention provides the ability to exploit spatial 
information to separate co-channel emitters for detection, classification, and 
identification. This is particularly applicable to detecting signals designed for a low 

1 0 probability of detection (LPD) or low probability of intercept (LPI) which may use 
ambient background electromagnetic radiation and known co-channel emitters as a 
means of concealment. The spatial fourth-order cumulant matrix-pencil based blind 
source separation technique in accordance with the present invention provides the 
capability to blindly separate unknown co-channel emitters that may be near or below 

1 5 the noise floor of the individual sensors. 

Although illustrated and described herein with reference to certain specific 
embodiments, the BSS technique as described herein is nevertheless not intended to 
be limited to the details shown. Rather, various modifications may be made in the 
details within the scope and range of equivalents of the claims and without departing 
20 from the spirit of the invention. 
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Claims 

What is claimed is: 

1 . A method for separating a plurality of signals provided by a respective 
plurality of sources and received by an array comprising a plurality of elements, said 

5 method comprising: 

generating a separation matrix as a function of: 

time differences between receipt of said plurality of signals by said 
plurality of elements; and 

a spatial fourth order cumulant matrix pencil; and 

10 multiplying said separation matrix by a time series matrix representation of 

said plurality of signals. 

2. A method in accordance with claim 1, wherein said spatial fourth order 
cumulant matrix pencil is a function of a spatial fourth order cumulant matrix. 

3. A method in accordance claim 2, wherein said spatial fourth order cumulant 
1 5 matrix is in accordance with the following equation: 

N 

Cj(r 1 ,r 2 ,r 3 ) = 5]Cti/w[x;(^r 1 )x l .(r-r 2 )x(r)x H (t - r 3 )] f where in: 

i=i 

C* (r, ,r 2 ,r 3 ) is said spatial fourth order cumulant matrix having a first 

time lag, xi, a second time lag, X2, and a third time lag, T3, each time lag 
being indicative of a time delay from one of said plurality of sources to 
20 one of said plurality of elements; 

M is indicative of a number of a number of elements in said array; 

Cum [xi*(t-xi)Xi(t-X2)x(t)x H (t-T3)] is a cumulant operator on arguments 
Xi\t.xi)xi(t-X2)x(t)x H (t-x 3 ); 

t is a variable representing time; 

25 Xj*(t-xi) represents a complex conjugate of one of said plurality of 

signals from an I th source at time t-xi; 

xj(t-X2) represents one of said plurality of signals from an i source at 
time t-xi ; 
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x(t) is a vector representation of said plurality of signals; and 

x H (t-i3) represents the Hermitian transpose of x(t-i3). 

4. A method in accordance with claim 1 , wherein: 

at least two of said plurality of elements have non-identical beam patterns. 

5 5. A method in accordance with claim 1 wherein: 

all of said plurality of elements have non-identical beam patterns. 

6. A method in accordance with claim 1 , wherein said step of generating said 
separation matrix comprises performing a generalized eigenvalue analysis of said 
spatial fourth order cumulant matrix pencil. 

10 7. A method in accordance with claim 1 , further comprising: 

estimating said spatial fourth order cumulant matrix pencil as a function of 
selected ones of said time differences; 

determining non-zero finite eigenvalues for said spatial fourth order cumulant 
matrix pencil; 

1 5 determining a number of said finite eigenvalues that are distinct; 

determining a multiplicity of each of said distinct finite eigenvalues; 

calculating linearly independent eigenvectors for each of said distinct finite 
eigenvalues; 

for each eigenvalue having a multiplicity equal to one, calculating a 
20 normalization factor and generating a respective separation vector as a 

function of said normalization factor and an eigenvector corresponding to said 
eigenvalue having a multiplicity equal to one; 

for each repeated eigenvalue, generating a linearly independent set of 
separation vectors using the multiplicity of eigenvectors belonging to each 
25 repeated eigenvalue; and 

generating said separation matrix as a function of said separation vectors. 

8. A method in accordance with claim 1 , wherein said time differences are not 
equal to zero. 
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9. A method in accordance with claim 1 , further comprising calculating an 
efficiency of separating a signal from said plurality of signals, wherein said efficiency 
is a function of a ratio of a power of a separated signal and a power of a signal from a 
respective source. 

5 10. A method in accordance with claim 9, wherein said efficiency is in accordance 
with the following equation: 

Qj=Sj/Pj, wherein: 

Cj is indicative of said separation power efficiency for a j source of 
said plurality of sources; 

10 Sj is indicative of a power of a separated signal from said source; 

and 

th 

Pj is indicative of a normalized power of a signal from said j source. 

11. A computer readable medium encoded with a computer program code for 
directing a processor to separating a plurality of signals provided by a respective 

1 5 plurality of sources and received by an array comprising a plurality of elements, said 
program code comprising: 

a first code segment for causing said processor to generate a separation matrix 
as a function of: 

time differences between receipt of said plurality of signals by said 
20 plurality of elements; and 

a spatial fourth order cumulant matrix pencil; and 

a second code segment for causing said processor to multiply said separation 
matrix by a time series matrix representation of said plurality of signals. 

12. A computer readable in accordance with claim 1 1 , wherein: 

25 said spatial fourth order cumulant matrix pencil is a function of a spatial fourth 

order cumulant matrix being a summation of steering vector outer products 
scaled by an individual source signal's fourth order cumulant; and 

said steering vector is indicative of respective phase delays between each 
element of said plurality of elements. 

30 
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13. A computer readable medium in accordance claim 12, wherein said spatial 
fourth order cumulant matrix is in accordance with the following equation: 

N 

Cj (r, ,t 2 , t 3 ) ^ ]T Cw/n [xj (f - r, ) x, (f - r 2 ) x (t ) x H (/ - r 3 )] 5 wherein: 

C A X (r, ,r 2 , r 3 ) is said spatial fourth order cumulant matrix having a first 

5 time lag, xj, a second time lag, x 2 , and a third time lag, X3, each time lag 

being indicative of a time delay from one of said plurality of sources to 
one of said plurality of elements; 

M is indicative of a number of a number of elements in said array; 

Cum [xi*(t-Ti)xi(t-T2)x(t)x H (t-t3)] is a cumulant operator on arguments 
1 0 x i \M,)x i (t-x 2 )x(t)x ll (t-T3); 

t is a variable representing time; 

Xi*(t-xi) represents a complex conjugate of one of said plurality of 
signals from an 1 th source at time t-xi; 

Xi(t-X2) represents one of said plurality of signals from an i source at 
15 timet-xi; 

x(t) is a vector representation of said plurality of signals; and 

x H (t-X3) represents the Hermitian transpose of x(t-X3). 

14. A computer readable medium in accordance with claim 1 1 , wherein: 

at least two of said plurality of elements have non-identical beam patterns. 

20 15. A computer readable medium in accordance with claim 1 1 , wherein: 

all of said plurality of elements have non-identical beam patterns. 



16. A computer readable medium in accordance with claim 1 1 , said program code 
further comprising: 

25 a third code segment for causing said processor to perform a generalized 

eigenvalue analysis of said spatial fourth order cumulant matrix pencil. 
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17. A computer readable medium in accordance with claim 1 1 , said program code 
further comprising: 

a fourth code segment for causing said processor to estimate said spatial fourth 
order cumulant matrix pencil as a function of selected ones of said time 
5 differences; 

a fifth code segment for causing said processor to determine non-zero finite 
eigenvalues for said spatial fourth order cumulant matrix pencil; 

a sixth code segment for causing said processor to determine a number of said 
finite eigenvalues that are distinct; 

10 a seventh code segment for causing said processor to determine a multiplicity 

of each of said distinct finite eigenvalues; 

an eighth code segment for causing said processor to calculate linearly 
independent eigenvectors for each of said distinct finite eigenvalues; 

for each eigenvalue having a multiplicity equal to one, a ninth code segment 
1 5 for causing said processor to calculate a normalization factor and generating a 

respective separation vector as a function of said normalization factor and an 
eigenvector corresponding to said eigenvalue having a multiplicity equal to 
one; 

for each repeated eigenvalue, a tenth code segment for causing said processor 
20 to generate a separation vector as a function of an eigenvector corresponding 

to said repeated eigenvalue; and 

an eleventh code segment for causing said processor to generate said 
separation matrix as a function of said separation vectors. 

18. A computer readable medium in accordance with claim 1 1, wherein said time 
25 differences are not equal to zero. 

19. A computer readable medium in accordance with claim 11, said program code 
further comprising: 

an efficiency calculation code segment for calculating an efficiency of 
separating a signal from said plurality of signals, wherein said efficiency is a 
30 function of a ratio of a power of a separated signal and a power of a signal 

from a respective source. 
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20. A computer readable medium in accordance with claim 1 9, wherein said 
efficiency is in accordance with the following equation: 

tj=Sj/Pj t wherein: 

Q is indicative of said separation power efficiency for a 3 th source of 
5 said plurality of sources; 

Sj is indicative of a power of a separated signal from said j source; 
and 

Pj is indicative of a normalized power of a signal from said j source. 

21 . A system for separating a plurality of signals provided by a respective 
10 plurality of sources, said system comprising: 

a receiver for receiving said plurality of signals and for providing received 
signals; and 

a signal processor for receiving said received signals, generating a separation 
matrix, and multiplying said separation matrix by a time series matrix 
1 5 representation of said received signals, wherein: 

said separation matrix is a function of time differences between receipt 
of said plurality of signals by said receiver and a function of a spatial 
fourth order cumulant matrix pencil. 

22. A system in accordance with claim 21, wherein said receiver comprises a 
20 plurality of elements configured to form an array. 

23. A system in accordance with claim 22, wherein: 

said spatial fourth order cumulant matrix pencil is a function of a spatial fourth 
order cumulant matrix being a summation of steering vector outer products 
scaled by an individual source signal's fourth order cumulant; and 

25 said steering vector is indicative of respective phase delays between each 

element of said plurality of elements. 
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24. A system in accordance claim 23, wherein said spatial fourth order cumulant 
matrix is in accordance with the following equation: 

N 

c£(T,,r 2 ,r 3 ) S ^ 

C\ (r, , t 2 , r 3 ) is said spatial fourth order cumulant matrix having a first 

5 time lag, Xj, a second time lag, x 2 , and a third time lag, 13, each time lag 

being indicative of a time delay from one of said plurality of sources to 
one of said plurality of elements; 

M is indicative of a number of a number of elements in said array; 

Cum [xi*(t-Ti)Xi(t-T2)x(t)x H (t-t3)] is a cumulant operator on arguments 
1 0 Xi*(t-x, )xi(M 2 )x(t)x H (t-T 3 ); 

t is a variable representing time; 

Xi*(t-Ti) represents a complex conjugate of one of said plurality of 
signals from an 1 th source at time t-X] ; 

Xi(t-x 2 ) represents one of said plurality of signals from an i* source at 
1 5 time t-xi ; 

x(t) is a vector representation of said plurality of signals; and 

x H (t-X3> represents the Hermitian transpose of x(t-X3). 

25. A system in accordance with claim 22, wherein: 

at least two of said plurality of elements have non-identical beam patterns. 

20 26. A system in accordance with claim 22, wherein: 

all of said plurality of elements have non-identical beam patterns. 

27. A system in accordance with claim 21, wherein said time differences are not 
equal. 
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28. A system in accordance with claim 21, said signal processor comprising: 

a matrix pencil estimation portion for estimating said spatial fourth order 
cumulant matrix pencil as a function of selected ones of said time differences; 

a non-zero finite eigenvalue determination portion for determining non-zero 
5 finite eigenvalues for said spatial fourth order cumulant matrix pencil; 

a number of distinct eigenvalue determination portion for determining a 
number of said finite eigenvalues that are distinct; 

a multiplicity determination portion for determining a multiplicity of each of 
said distinct finite eigenvalues; 

10 a linearly independent eigenvector calculation portion for calculating linearly 

independent eigenvectors for each of said distinct finite eigenvalues; 

a normalization factor portion for calculating, for each eigenvalue having a 
multiplicity equal to one, a normalization factor and generating a respective 
separation vector as a function of said normalization factor and an eigenvector 
1 5 corresponding to said eigenvalue having a multiplicity equal to one; 

a separation vector generation portion for generating, for each repeated 
eigenvalue, a separation vector as a function of an eigenvector corresponding 
to said repeated eigenvalue; and 

a separation matrix generating portion for generating said separation matrix as 
20 a function of said separation vectors. 

29. A system in accordance with claim 21 further comprising a separation power 
efficiency portion for calculating an efficiency of separating said plurality of signals 
in accordance with the following equation: 

Q=Sj/Pj ? wherein: 

25 Q is indicative of said separation power efficiency for a source of said plurality of 
sources; 

Sj is indicative of a power of a separated signal from said j source; and 
Pj is indicative of a normalized power of a signal from said j source. 
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